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SECTION  1 


INTRODUCTION 

The  proper  simulation  of  plasma  transport  in  three 
dimensions  is  a  crucial  issue  in  the  construction  of  a 
numerical  model  for  the  high  altitude  nuclear  environment 
(HANE) .  As  a  rule,  the  digital  simulation  of  the  system 
of  spatial  equations,  which  describe  the  physics,  requires 
that  all  fluctuations  of  the  simulated  quantities  be 
resolved.  In  general,  HANE  simulations  cannot  meet  this 
requirement  because  the  range  of  relevant  scale  sizes  can 
extend  to  six  orders  of  magnitude.  The  practical  limita¬ 
tions  of  computer  storage  and  cost  will,  typically,  con¬ 
strain  a  three  dimensional  numerical  code  to  a  scale  size 
range  of  one  to  one-and-one  half  orders  of  magnitude  in 
each  dimension.  At  this  sampling  rate,  only  the  larger 
structures  can  be  resolved.  Historically,  HANE  plasma 
simulation  algorithms,  as  employed  in  codes  such  as  MELT, 
spatially  resolve  the  large  scales,  but  are  forced  to  ignore 
a  vast  range  of  small  scales.  It  is,  of  course,  generally 
acknowledged  that  the  small  scales  are  important  to  satellite 
and  IR  systems  effects  modelling,  but  their  interaction  with 
the  large  scale  flow  has  been  generally  neglected  in  codes 
such  as  MELT.  The  importance  of  small  scale  structure 
effects  has  motivated  the  development  of  a  theory  of  micro¬ 
structure  generation  and  transport. 

In  recent  years,  a  theory  of  microstructure  generation 
and  transport,  based  on  the  gradient  drift  mechanism,  has 
been  developed  which  is  capable  of  addressing  the  scale  size 
range  problem  in  computer  simulation.  The  numerical 
algorithms  which  derive  from  this  microstructure  theory  have 
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been  incorporated  into  the  nuclear  effects  simulation 
SCENARIO  code.  In  terms  of  general  code  architecture, 
SCENARIO  may  be  viewed  as  an  engineering  version  of  MELT. 
The  principal  purpose  of  the  code  is  to  provide  three- 
dimensional  plasma  environments  for  subsequent  systems 
analysis  of  satellite  communications  effects.  In  the 
present  context,  the  code  is  serving  as  a  research  vehicle 
for  developing  and  applying  the  microstructure  theory. 

This  paper  focuses  upon  the  cross fie Id  transport  imple¬ 
mented  in  the  code  which  represents  the  recent  advances 
in  the  theory. 

Microstructure  theory  can  be  regarded  as  the  applica¬ 
tion  of  stochastic  techniques  to  spatial  simulation 
problems.  The  application  of  stochastic  theory  to  the 
HANE  plasma  problem  originated  with  research  by  Workman 
and  Chu  in  1975.  Utilizing  the  principle  of  conservation 
of  plasma  power  (the  square  of  electron  density) ,  the  early 
work  tied  the  physical  space  numerical  concepts  into  the 
Fourier  space  data  base  (as  expressed  in  a  power  spectral 
density  curve) .  The  first  theoretical  and  numerical  model 
concentrated  on  the  evolution  of  plasma  scale  size  in 
Fourier  space  (Workman  and  Chu,  1978) .  This  paper  has 
stimulated  other  interesting  research  into  the  Fourier 
space  evolution  of  plasma  power,  notably  the  outstanding 
experimental  work  by  Pongratz  and  Fitzgerald  (1981) . 

The  early  microstructure  work  was  followed  by  the 
development  of  transport  algorithms  for  the  evolution  of 
both  mean  plasma  density  and  density  variance.  This 
modelling  accounted  for  the  effect  of  microstructure 
generation  on  gross  plasma  structures  and  flow  (Workman, 
Chu,  and  Ferrante,  1979).  The  important  concept  of  the 
flow  of  plasma  power  from  macrostructure  scales  to 
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microstructure  scales  in  a  conserved  manner  was  extensively 
developed  in  this  work.  The  introduction  of  the  "Haerendel 
velocity"  approximation  as  a  suitable  vehicle  for  furthering 
the  theoretical  development  was  made  by  Workman  in  connection 
with  the  work. 

A  significant  advance  in  the  analytic  development  of 
microstructure  theory  towards  a  more  formal  and  rigorous 
stochastic  approach  was  undertaken  by  Stagat  in  1979.  The 
"smooth  current"  approximation  in  his  work  provided  both  a 
more  complete  statement  of  the  Haerendel  velocity  model 
and,  more  important,  laid  the  groundwork  for  a  self  consistent 
and  elegant  mathematical  formulation  of  the  general  problem. 
Further  developments  {as  reported  in  Stagat,  et  al.,  1982) 
have  formed  the  starting  point  for  all  subsequent  work. 

This  paper  is  directed  at  recent  developments  in 
microstructure  theory,  but  is  organized,  in  its  derivation 
of  results  from  basic  principles,  to  be  a  self  contained 
document.  That  is,  while  the  reader  is  encouraged  to  study 
the  earlier  work  for  the  physical  insight  that  has  led  to 
current  theory,  the  objective  here  is  to  present  a  complete 
development  of  the  work  as  it  now  stands. 

The  application  of  stochastic  theory  to  the  equations 
for  the  evolution  of  2-dimensional  plasmas,  as  discussed  in 
Section  2,  leads  to  further  understanding  of  the  implications 
of  the  "Haerendel  velocity"  and  the  "smooth  current"  or 
"smooth  J"  assumptions  in  terms  of  the  correlation  properties 
assumed  between  the  stochastic  quantities  of  the  problem 
and  in  terms  of  the  modifications  to  the  effective  conduc¬ 
tivity  that  are  assumed  to  be  caused  by  microstructure. 
Two-dimensional  algorithms  are  developed  based  on  this 
theory  and  are  applied  to  the  2-dimensional  evolution  both 
of  a  barium  ion  cloud  and  a  nuclear  plasma  under  the 
influence  of  radial  heave  winds. 
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In  Section  3  the  rigorous  approach  used  to  develop 
the  2-diroensional  algorithms  in  Section  2  are  extended  to 
the  problem  of  3-dimensional  transverse  transport.  The 
transverse  transport  aspects  of  3-dimensional  late-time 
HANE  plasmas  are  very  similar  to  the  transport  and  evolu¬ 
tion  of  2-dimensional  plasmas.  Assumptions  about  the 
structure  correlation  properties  in  the  parallel-to-the- 
field  direction  are  needed  in  the  3-dimensional  HANE 
problem  to  remove  this  dimension  from  the  problem,  but 
once  this  removal  has  been  accomplished  the  two  problems 
are  nearly  identical.  Different  transport  algorithms 
arise  from  different  assumptions  and  several  possible 
assumptions  are  considered  in  Section  3.  Algorithms 
based  on  one  promising  set  of  assumptions  have  been 
developed  and  implemented  in  a  version  of  SCENARIO. 
Results  from  a  trial  run  are  presented. 

Section  4  presents  conclusions  as  well  as  other 
issues  into  which  further  extension  of  the  theory  is 


SECTION  2 


MICROSTRUCTURE  THEORY  APPLIED  TO  2-DIMENSIONAL 
PLASMA  EVOLUTION 

2-1  CLASSICAL  SPLIT  STEP  ALGORITHM  FOR  2-DIMENSIONAL 

PLASMA  EVOLUTION 

Classical  modelling  of  2-dimensional  plasma  evolution 
is  performed  using  split-step  algorithms  that  are  characterized 
by  the  computational  flow  schematized  in  Figure  2-1.  Two  key 
steps  seen  in  the  figure  are  the  potential  solving  step  and 
the  convection  step.  Given  a  specification  for  the  spatially 
varying  conductivity,  the  potential  step  finds  the  current 
flow  and  electric  field  patterns  produced  by  plasma  driving 
forces.  The  plasma  motion  is  described  by  electron  Hall 
drift  (E  x  B  motion)  and  thus  electric  field  determination 
is  equivalent  to  plasma  velocity  determination.  Given  the 
results  of  the  potential  equation  solution,  the  plasma  flow 
is  known  and  can  be  used  to  advance  plasma  density  parameters 
through  an  increment  in  time  in  the  convection  step.  The 
changes  in  plasma  density  cause  corresponding  changes  in  the 
plasma  conductivity  which  are  provided  to  the  potential  solver 
through  the  conductivity  relation  for  the  start  of  the  next 
split-step  cycle. 

For  the  classical'  case  the  conductivity  relation  is  a 
simple  proportionality  and  is,  thus,  not  often  singled  out  in 
classical  split-step  descriptions  as  it  is  in  Figure  2-1. 
Portraying  it  as  an  explicit  step  in  this  instance  provides 
an  early  insight  into  the  importance  this  relation  has  in 
microstructure  algorithms.  Microstructure  algorithms  for 
plasma  evolution  can  be  regarded  as  an  extension  of  the 
split-step  scheme  shown  in  Figure  2-1,  as  will  be  discussed 
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Potential  Eq:  v*Cv<t>«  v-  C(V  xiz) 
Convection  Eq:  &tN  +v*NU  ■  0 
Conductivity  Eq:  C  ■  eN/hB 


Figure  2-1.  Classical  2-d  split  step  algorithm 
for  plasma  evolution. 
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later.  It  is  appropriate  to  discuss  in  this  section  the 
background  of  the  relations  behind  the  steps  portrayed. 

A  key  physical  parameter  in  plasma  evolution  behavior 
is  the  plasma  current  density,  J.  Plasma  currents  are 
driven  by  imbalances  between  the  forces  acting  on  the 
electron  fluid  and  the  forces  acting  on  the  ion  fluid. 

Currents  are  driven  not  only  by  the  conventional  electro¬ 
magnetic  field  forces  but  also  by  neutrdl  wind,  pressure, 
gravity,  and  centrifugal  acceleration  forces.  For  the  late 
time  HANE  plasma  situation  the  electromagnetic  field  forces 
and  the  neutral  wind  forces  frequently  dominate  plasma 
behavior.  For  the  sake  of  simplicity  this  report  will 
concentrate  on  these  two  forces.  The  results  obtained,  it 
should  be  noted,  are  easily  generalizable  to  a  more  complete 
list  of  HANE  forces. 

The  ionospheric  current  response  to  a  combination  of 
electric  and  neutral  wind  forces  is  approximated  in  this 
report  with  the  relation 

J-fi(l+Vx  **)  (2-1) 

where  e  is  the  electronic  charge,  N  is  the  plasma  density, 
y  is  the  ionic  mobility  through  the  neutral  air  background, 

V  is  the  neutral  wind,  E  is  the  electric  field  and  B  is 
the  magnetic  field  intensity  and  where  the  magnetic  field  is 
assumed  to  be  constant  and  directed  in  the  z-direction.  For 
the  two-dimensional  evolution  problem  the  quantities  N,  E 
and  V  are  assumed  to  vary  only  in  the  x-  and  y-directions. 
The  quantity  i  represents  the  unit  normal  vector  in  the 

Z 

z-direction.  Note  that  eN/yB  can  be  identified  as  the 
plasma  conductivity. 
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The  above  current  response  is  found  from  an  evaluation 
of  the  difference  between  ion  and  electron  fluxes  in  response 
to  electromagnetic  field  and  wind  forces.  The  electron 
velocity  is  assumed  to  be  essentially  independent  of  the  wind 
forces.  The  electron  mobility  in  the  background  air  is  also 
assumed  to  be  very  high  implying 

U  -  E  x  i  /B  (2.2) 

e  z 

where  U  is  the  electron  velocity. 

The  ion  velocity  is  evaluated  from  the  momentum  or 
force  balance  equation.  A  force  balance  equation  between  the 
ion  drag  in  a  neutral  gas  and  an  imposed  electric  field 
produces  the  simple  proportionality  relation 


U±  *  pE  (2.3) 

where  IL  is  the  ion  velocity  and  where  again  y  is  the 
ionic  mobility  in  the  gas.  If  the  neutral  gas  is  not  at  rest 
but  is  moving  with  velocity  V  then  the  drag  relation 
becomes 

Ui  -  V  »  yE  (2.4) 

Heuristically  it  is  not  expected  that  ions  could  differentiate 
between  electromagnetic  forces  that  are  purely  electrostatic 
and  forces  that  include  magnetic  fields.  This  argument 
justifies  use  of  the  above  expression  for  the  more  general 
case  of  having  a  magnetic  field  present  with  the  electric 
field  E  in  (2.4)  substituted  with  E  +  (7^  *  B  giving 

0i  -  V  »  y (E  +  Ui  x  B)  (2.5) 
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Algebraic  manipulation  of  this  equation  to  isolate 
produces  the  well  known  result 

n  =  (UB)2  E  .  UB  E 

1  1  +  (UB)2  B  Z  1  +  (UB)2  B 


UB 


1  +  (UB) 


T  V  X  i2  + 


1  +  (UB) 


(2.6) 


2 

The  ion  velocity  for  large  uB  is  essentially  E  x  B/B 
and  thus  very  close  to  the  electron  velocity.  The  small 
differences  between  the  ion  velocity  and  the  electron  velocity 
are  the  source  of  plasma  currents.  The  difference  between  ion 
and  electron  velocity  can  be  seen  to  be 


U.  -  U  = 
i  e 


1  +  (UB) 


E 

7  B  X 


+ 


1  +  (UB) 


E 

7  B 


(2.7) 

+  - -  v  x  i  +  - - - *•  V 

1  +  (UB) *  Z  1  +  (UB)J 

which  for  large  yB  produces  the  approximate  current  relation 
presented  earlier  in  the  discussion. 

As  mentioned  previously  plasma  motion  is  taken  as 
synonymous  to  electron  motion  which  is  affected  only  by  the 
electric  field.  The  mechanism  by  which  wind  forces  move  a 
plasma  irregularity  must  thus  be  through  the  electric  fields 
it  generates  as  there  are  no  explicit  wind  terms  in  the 
electron  velocity.  The  mechanism  by  which  the  wind  generates 
electric  fields  is  through  the  current  that  it  drives.  Where 
the  wind  current  encounters  gradients  in  density  its  flow  is 
modified  so  that  it  is  no  longer  divergenceless.  Charge 
density  source  regions  and  sink  regions  are  created  which 
build  up  charge.  This  build  up  continues  over  the  period  of 
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typically  milliseconds  until  the  associated  electric  fields 
build  up  to  a  point  where  electrically  driven  currents  bleed 
off  the  charge  sources  and  sinks.  The  balance  between  the  wind 
current  produced  charge  sources  and  sinks  and  the  electric 
field  driven  current  is  specified  by  the  relation 

V  •  C  E/B  »  -  V  •  CV  X  i  (2.8) 

z 

where  C,  the  conductivity,  is  given  by  eN/pB.  Thus  the 
plasma  flow  in  response  to  applied  wind  forces  is  established 
through  the  creation  of  a  polarization  field  which  can  be 
found  through  the  above  equation.  The  additional  assumption 
that  V  *  E  =  -3B/3t  *  0  is  made  for  late  time  electrostatic 
flow  which  facilitates  solution  for  E  in  terms  of  a 
potential  function ,  4> .  The  resulting  relations  are 


E/B  *  -  7* 

V  •  CV<J>  »  V  •  C  V  X  iz 


(2.9) 


where  C  and  V  are  known  functions  of  x  and  y  and 
where  <J>  is  to  be  found.  The  potential  <P  is  found  through 
these  relations  in  the  potential  step  algorithm. 

The  convection  step  algorithm  uses  the  plasma 
continuity  equation  with  the  known  plasma  flow  to  update 
the  plasma  density.  The  continuity  relation  is  given  by 


+  v  •  NU  «  0  (2.10) 

where  U  is  the  plasma  flow  velocity  given  by  Ex  iz/B. 

The  continuity  equation  is  basically  a  restatement  of  the 
conservation  of  electrons.  The  only  source  of  local  density 
changes  is  assumed  to  be  variations  of  the  fluxes  into  and 


out  of  a  region.  The  continuity  expression  can  also  be 
written  in  terms  of  the  current  density  J 


3N 

at 


+ 


i  +  N  V 


=  0 


(2.11) 


This  form  of  the  continuity  expression  is  sometimes  more 
useful  for  microstructure  transport  considerations. 

The  background  relations  for  a  complete  cycle  of  the 
split  step  algorithm  have  thus  been  presented.  The  relation 
between  the  local  density  and  the  plasma  conductivity  is 
obtained  from  the  current  relationship  with  applied  electric 
field  and  wind  forces.  The  plasma  motion  determined  by 
the  potential  step  is  found  by  balancing  the  plasma  charge 
buildup  and  depletion  caused  by  wind  current  divergences 
with  electric  current  flows  produced  by  electric  fields. 

The  flow  quantities  calculated  are  then  used  with  the  conserva¬ 
tion  principle  to  update  plasma  densities.  The  next  section 
discusses  how  these  steps  are  modified  to  account  for  the 
presence  of  microstructure. 


2-2  APPLICATION  OF  STOCHASTIC  THEORY  TO  THE  TWO-DIMENSIONAL 

SPLIT-STEP  ALGORITHM 

In  the  previous  section  a  description  of  the  classical 
split-step  approach  to  the  simulation  of  plasma  evolution  is 
presented.  In  this  section  the  basic  refinements  provided 
by  the  application  of  microstructure  theory  to  the  plasma 
transport  algorithms  are  discussed. 

The  application  of  microstructure  theory  to  the  plasma 
evolution  problem  is  motivated  by  one  primary  concern,  namely 
fine  scale  plasma  structure.  Fine  scale  plasma  structure 
with  sizes  down  to  the  order  of  100  meters  is  produced  over 
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spatial  extents  of  hundreds  of  kilometers  by  high  altitude 
nuclear  explosions.  Its  location  intensity  and  spatial 
characteristics  are  important  for  the  prediction  of  satellite 
communications  effects.  Because  of  the  complex  and  seemingly 
chaotic  nature  of  the  HANE  structure,  stochastic  descriptions 
of  it  are  used  to  predict  satellite  communications  effects. 

Only  a  few  statistics  of  the  fine  scale  structure  are 
required  to  provide  reliable  estimates  of  satcom  system  per¬ 
formance  and  these  statistics  are  believed  to  be  nearly  con¬ 
stant  over  many  kilometer  regions. 

It  is  virtually  impossible  to  grid  HANE  numerical 
computations  fine  enough  to  resolve  the  fine-scale  structure 
details  because  of  storage  limitations.  Classical  techniques 
merely  grid  as  fine  as  costs  and  storage  limitations  permit 
and  ignore  the  presence  of  plasma  density  and  velocity 
fluctuations  varying  faster  than  the  spatial  Nyquist  rate. 

These  approaches  do  not  specify  where  the  very  important  fine 
scale  density  fluctuations  exist  and  they  ignore  any  possible 
effects  of  the  fine  scale  structure  on  transport.  Micro¬ 
structure  theory  was  developed  to  address  these  deficiencies 
in  the  classical  techniques. 

Microstructure  theory  treats  the  plasma  evolution 
equations  on  a  stochastic  basis.  Numerical  algorithms  are 
used  with  sampling  grids  that  have  roughly  the  same  resolution 
as  the  classical  problem  and  that  cannot  resolve  the  fine 
scale  structure  on  a  deterministic  basis.  Rather  than  ignore 
the  presence  of  the  fine  scale  fluctuations,  the  microstructure 
approach  treats  structure  unresolved  by  the  sampling  grid  on 
a  stochastic  basis  and  describes  its  characteristics  through 
its  statistics.  (Structure  unresolved  by  the  sampling  grid 
when  approached  on  a  stochastic  basis  is  what  is  meant  by  the 
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term  ’micro structure. ' )  By  adopting  a  stochastic  approach 
to  the  problem,  transport  algorithms  for  structure  location 
and  intensity  can  be  developed  for  satcom  effects  predictions; 
the  influence  of  fine  scale  structure  upon  bulk  plasma  motion 
can  be  investigated;  and  an  understanding  of  the  implicit 
assumptions  behind  the  classical  approach  can  be  developed. 

Microstructure  theory  provides  refinements  to  the 
classical  plasma  evolution  algorithms  which  are  evident  in 
a  comparison  of  Figure  2-1  with  Figure  2-2,  the  latter  of 
which  illustrates  the  microstructure  split-step  algorithms. 

The  potential  solving  steps  of  both  algorithms  are  the  same. 

Both  produce  plasma  flow  parameters  for  the  convection  step 
from  the  specified  spatial  dependence  of  conductivity.  The 
fundamental  difference  between  the  two  algorithms  is  that  the 
conductivity  specified  for  the  microstructure  algorithm  is  a 
function  of  multiple  statistical  parameters  instead  of  the 
one  parameter  of  plasma  density  used  in  the  classical  algorithm. 
Further,  in  order  to  supply  the  spatial  dependence  of  these 
multiple  parameters  it  is  necessary  for  the  microstructure 
convection  algorithm  to  convect  multiple  statistical  quantities. 

The  refinements  to  the  classical  algorithm  provided  by 
the  microstructure  approach  are  based  upon  the  premises  that 
fine  scale  structure  changes  the  conductivity  from  the  value 
determined  by  the  mean  plasma  density  in  a  cell  and  that  by 
specifying  statistical  quantities  in  addition  to  the  mean  a 
better  approximation  to  the  actual  conductivity  may  be  obtained. 
Specifying  one  or  two  additional  quantities  requires  a  loss  in 
the  linear  scale  resolution  which  is  small  in  comparison  to 
the  gain  in  the  accuracy  of  bulk  plasma  convection  and  which  is 
virtually  inconsequental  in  view  of  the  distinct  advantage  of 
being  able  to  predict  structure  location  and  intensity. 


Conductivity 
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Algorithms 


<  J  >  and  V 

OR 
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Potential  Eq:  v-Cv*  =  v*C (V  xiz> 

Convection  Eq:  dT<F(N)>  +v.  <f(N)U>  =  0 

The  term  <NU>  is  obtained  from<J>. 

Determining  <f(n)U>  requires  assumptions 
about  correlation  between  <t> ,  J  and  N. 

Conductivity  given  in  terms  of  Statistics  of  N 


Figure  2-2.  Microstructure  2-d  split  step  algorithm 
for  plasma  evolution. 
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In  the  microstructure  approach  the  basic  quantities 
of  plasma  density,  plasma  velocity,  electrostatic  potential 
and  plasma  current  are  treated  as  stochastic  processes  which 
are  interrelated  and  advanced  through  time  with  the  classical 
differential  equations  understood  in  a  stochastic  sense. 

The  driving  force  of  the  neutral  wind  is  taken  to  be 
deterministic.  The  stochastic  quantities  are  assumed  to  have 
homogeneous  statistics  over  each  cell.  The  statistics  are 
treated  as  inhomogeneous  over  the  simulation  dimension  and 
it  is  the  spatial  dependence  of  the  statistics  that  the 
approach  tries  to  predict. 

In  the  potential  equation  and  in  the  convection  equations 
for  the  statistics  the  stochastic  quantities  are  strongly 
coupled  together.  As  a  result  the  nature  of  the  correlation 
between  them  is  important  for  characterizing  rheir  evolution. 
Currently  the  correlation  properties  between  plasma  quantities 
are  not  predicted  quantities.  They  are  determined  instead 
through  basic  assumption  (justified  perhaps  through  physical 
intuition  and  to  some  extent  by  numerical  investigation) . 

It  is  perhaps  surprising  that  the  classical  algorithm, 
which  simply  ignores  the  presence  of  fine-scale  fluctuations, 
can  be  justified  from  a  microstructure  perspective  by  making 
the  simple  assumption  that  the  electric  field  and  the  density 
fluctuations  are  uncorrelated.  This  assumption  allows  the 
conductivity  relation  for  the  mean  current  <J>  to  be 
derived  straightforwardly  by  taking  the  expectation  (denoted 
by  angle  irons  <  >)  of  both  sides  of  the  current  relation 
(Eq.  2.1)  presented  in  Section  2.1: 


<J>  *  1  <N  |>  +  4  <N>  V  x  i 
uB  B  uB  z 


e<N>  /< E> 
yB  \  B 


+  V  x 


(2.12) 
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The  latter  form  of  this  equation  has  incorporated  the  assumed 
uncorrelated  nature  of  the  fluctuations  of  the  density  N  and 
the  electric  field  E. 

The  mean  plasma  density  <N>  is  seen  to  play  the  same 
role  as  the  density  N  does  in  the  classical  conductivity 
relation.  The  convection  relation  for  the  mean  plasma  density 
<N>  can  also  be  found  to  be  the  same  as  that  for  the  classical 
density  by  taking  the  expectation  of  the  continuity  equation 


3<N>  _ 


<V 


N  E>  x  i  /B 
z 


a  7  •  <N>  <E>  X  i  /B 

z' 


(2.13) 


Thus,  microstructure  theory  can  justify  the  classical  con¬ 
vection  algorithms  under  the  assumption  that  the  density  and 
the  electric  field  are  uncorrelated  and  by  making  a  corre¬ 
spondence  between  the  classically  specified  density  and  the 
mean  density. 

It  may  be  of  further  interest  to  note  that  microstructure 

theory  techniques  can  be  applied  to  predict  the  convection  of 

the  variance  of  the  plasma  fluctuations  using  a  similar 

2 

assumption.  If  both  N  and  N  are  assumed  to  be 

uncorrelated  with  the  electric  field  then  it  is  found  that 
2 

<N  >  is  convected  with  the  same  velocity  as  <N>  and 

2 

consequently  the  plasma  variance  E  moves  with  the  same 

velocity  as  <N>.  To  illustrate  the  convection  velocity  of 
2 

<N  >  note 

-  2N  |S  *  -2N  E  x  i  •  VN/B 

o  t  0U  Z 

*  -V  .  N2  E  X  i  /B 
z 

and ,  thus , 
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(2.14) 


3N2  _  S<N2>  ... 

<TT  'at 


“<7 


E  x  i  /B> 
z 


2 

=  -7  •  <N  >  <E>  x  i  /B 

7 

where  the  last  expression  incorporates  the  assumed  uncorrelated 

2  9 

nature  between  N  and  E.  The  convection  velocity  of  <N  > 
is  thus  <E>  x  i_/B  which  is  the  same  as  the  convection 

Z 

velocity  of  <N>. 

Thus  microstructure  theory  can  be  applied  to  the  plasma 
evolution  problem  with  assumptions  that  result  in  not  only 
the  convection  of  density  in  a  classical-like  manner  but 
also  the  convection  of  structure.  The  structure  convection 
under  these  assumptions  is  somewhat  trivial  in  that  it  moves 
with  the  same  velocity  as  the  mean  plasma  density  but  the 
exercise  is  worthwhile  for  understanding  the  relationships 
between  classical  and  other  microstructure  convection 
techniques . 

Alternative  assumptions  about  the  correlation  between 
stochastic  quantities  lead  to  alternative  conductivity 
relations.  Early  in  the  development  of  the  microstructure 
approach  for  SCENARIO,  Workman  utilized  the  "Haerendel 
velocity"  assumption  for  plasma  transport.  Stagat  developed 
identical  two-dimensional  transport  relations  utilizing  the 
"smooth  current  density"  approximation.  Couched  in  terms 
of  microstructure  theory  these  transport  relations  are  the 
immediate  consequence  of  assuming  that  the  inverse  of  plasma 
density  is  uncorrelated  with  the  current  density.  The 
conductivity  relation  is  thus  found  by 


I  *  A  (I  ♦ v  *  0 

and 

<i>  -  <J><1/H>  .  ^  .  V  X  i2] 

giving  , 

<J>  =  A  <1/N>-1  Hr  +  V  x  il  (2.15) 
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For  the  SCENARIO  algorithms  the  statistic  1/<1/N>  replaces 
the  mean  plasma  density  <N>  used  in  the  classical  conduc¬ 
tivity  algorithm.  In  general  this  statistic  can  be  shown 
to  be  less  than  or  equal  to  the  mean  density  with  larger 
differences  reflecting  the  presence  of  more  intense  structure. 
Note  that  this  conductivity  relation  thereby  reflects  an 
implicit  assumption  that  microstructure  modifies  conductivity. 
This  conductivity  relation  and  others  similar  to  it  will  be 
discussed  in  further  detail  in  the  next  section. 

The  convection  algorithms  for  the  mean  density  is 
straightforwardly  obtained  by  taking  the  mean  of  the 
continuity  equation  in  its  current  density  form  as  presented 
in  the  previous  section/  namely  Eq.  (2.11),  i.e., 

||  .  -7  •(?  J  *  i,  ♦  «■  v) 
yielding  in  the  mean 

=  "V  *0?  <J>  X  iz  +  <N>  V)  (2.16) 

This  expression  for  the  convection  of  mean  density  is  applic¬ 
able  regardless  of  the  correlation  assumption  because  of  the 
absence  of  product  terms  of  stochastic  quantities  in  the 
stochastic  continuity  equation  as  given  in  the  above  form. 
Differences  in  the  mean  density  transport  and  evolution  can 
thus  be  attributed  directly  to  differences  in  the  mean 
current  density  pattern.  Current  density  patterns,  of  course, 
reflect  differences  in  conductivity  dependences  produced  by 
different  correlation  assumptions. 

An  effective  plasma  transport  velocity  can  be  defined 
such  that  the  convection  equation  for  the  mean  density  can 
be  written  as  follows: 


For  the  SCENARIO  approach  the  effective  velocity  can  be 
written  in  terms  of  the  electric  field  as 


0  *  U  +  F  W 

eff  ave  c 


(2.18) 


where  U  is  the  average  plasma  velocity  defined  by 

clV© 


Uave  =  <E>  x  VB  =  "<V^  x  iz  ; 


(2.19) 


where  W  is  the  average  ion-neutral  slip  velocity  defined 
by 

»  *  v  -  Dave  !  (2-20) 


and  where  F  is  a  normalized  structure  intensity  parameter 

C 

varying  from  zero  in  structureless  regions  to  near  one  in 
highly  structured  regions  defined  by 


=  1  - 


<N>  <: 


(2.21) 


The  continuity  equation  when  written  in  this  format  shows 
the  effective  transport  velocity  for  the  mean  plasma  to  be 
between  the  average  plasma  velocity  as  determined  from  the 
potential  equation  and  the  neutral  wind  velocity.  The 
equation  indicates  that  the  more  highly  structured  regions 
tend  to  move  more  with  the  neutral  wind  than  with  the  average 
plasma  velocity.  This  observation  has  been  made  in  past 
works  by  Workman  and  Stagat.  When  comparison  is  made  to 
classical  transport  velocities  it  should  be  kept  in  mind 
that  this  tendency  is  partially  offset  by  the  fact  that  the 
highly  structured  regions  will  be  regions  of  lower  conductivity 
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and  U  •  will  consequently  be  less  coupled  to  the  neutral 
ave 

wind.  When  investigating  the  differences  in  mean  plasma 

transport  for  alternative  conductivity  models  it  may  thus  be 

advisable  to  use  the  current  density  form  of  the  continuity 

equation  and  to  compare  current  densities. 

The  motion  and  growth  of  structure  in  the  SCENARIO 

approach  can  be  described  through  the  equation  for  the 

2 

plasma  variance  l  .  This  equation  is  derived  from 

2 

equations  describing  the  convection  of  <N  >  and  <N>  as 
follows: 

ai2  _  3(<N2>  -  <N>2) 

at  at 

-  -  2<H> 

a  -7  .  <N2  U>  +  2<N>  7  •  <N  U>  (2.22) 

2 

The  mean  flux  of  N  can  be  found  by  multiplying  the  flux 

expression  for  N  in  terms  of  the  current  by  N  and  by  then 

taking  the  mean.  The  flux  expression  for  N  and  the  subse- 

2 

quent  operations  to  find  <N  U>  are: 

N  U  =  El*  J  x  i  +  NV 

e  z 

N2  D  »  ^  NJ  x  i  +  N2  V 
e  z 

2  UB  2 

<N  U>  *  —  <NJ>  X  i  +  <N  >  V 
e  z 

*  <N><J>  x  i  +  <N2>  V  (2.23) 

G  z 

In  the  ‘latter  expression  another  correlation  assumption  is 

made  in  the  SCENARIO  approach  that  the  plasma  density  is  also 

2 

uncorrelated  with  the  current  density.  For  Z  the 
resulting  expression  is 
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id  =  .7 .  fk 

at  v  \j 


^  <N>  <J>  x  i  +  < 
e  z 


+  2<N>  V  *1^  <J>  X  i  +  <N>  V 
|^e  z 


N2>  V^j 
} 


(2.24) 


An  expression  for  the  convection  of  z*  in  terms  of  only 
2 

Z  and  <N>  is  desirable  for  use  in  numerical  algorithms. 

2 

With  some  algebraic  manipulation  the  terms  involving  <N  > 
can  be  removed  from  the  above  equation  yielding 


1! - v 

at 


<E2  V  -  <N>2  0eIf>  -  Ueff  •  7  <H>2  (2 _25| 


where  ueff  ^-s  again  the  effect  transport  velocity  of  the 
mean  plasma  density  as  defined  above  in  equation  (2.18). 

In  the  above  format,  the  equation  is  physically 
enlightening.  The  equation  describes  the  variance  as  a 
compressible  fluid  with  an  effective  motion  given  by 

v  ‘  <J?~  °8ff 

and  with  a  source  term  given  by 

2 


U 


eff 


V  <N> 


The  source  term  demonstrates  a  behavior  that  is  related  to 
the  gradient  drift  growth  of  structure. 

2 

Unlike  the  convection  equations  for  <N>  and  <N  > 

equation  (2.25)  cannot  be  put  into  a  form  where  the  entire 

right-hand  side  is  a  divergence  operation  on  a  flux. 

This  feature  of  the  equation  reflects  the  physical  behavior 

2 

of  gradient  drift  growth  of  structure.  The  total  Z  in 

a  numerical  simulation  need  not  remain  constant  because  of 

this  growth  term.  Howevej:  it  should  be  noted  that  numerical 

2 

algorithms  which  convect  Z  with  the  above  equation  and 
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<N>  with  its  convection  relation  do  come  close  to  conserving 
2 

the  total  <N  >  within  the  simulation.  This  fact  is 
important  when  considering  which  systems  of  statistics  should 

be  used  for  microstructure  algorithms,  e.g.,  choice  of  the 

2  2 
system  of  <N>  and  £  versus  the  system  of  <N>  and  <N  > 

Further  numerical  considerations  are  discussed  in  Section  2-4 

Thus  it  has  been  demonstrated  that  the  plasma  evolution 

equations  can  be  considered  on  a  stochastic  basis  in  order 

to  accommodate  the  presence  of  fine-scale  structure. 

Algorithms  can  be  derived  which  predict  the  location  of 

micro structure  and  its  influence  upon  plasma  evolution. 

These  algorithms  are  dependent  upon  the  assumptions  made 

about  which  quantities  are  uncorrelated.  The  classical 

plasma  evolution  algorithms  set  forth  in  the  previous  section 

can  be  justified  from  a  stochastic  standpoint  if  the  electric 

field  fluctuations,  and  hence  the  velocity  field  fluctuations 

are  assumed  uncorrelated  with  the  plasma  density  fluctuations 

It  is  observed  earlier  in  this  section  that  the  classical 

approach  can  be  extended  to  predict  the  location  of  structure 

Using  the  extended  classical  approach  structure  is  found  to 

be  frozen  in.  For  SCENARIO  an  alternative  set  of  stochastic 

variables  are  assumed  to  be  uncorrelated,  specifically,  the 

density  and  its  inverse  both  are  assumed  to  be  uncorrelated 

with  the  current  density  fluctuations.  Using  the  SCENARIO 

approach  plasma  structure  is  found  to  move  with  a  slightly 

different  velocity  than  the  mean  density  and  to  grow  with 

movements  through  density  gradients.  Justification  for  the 

SCENARIO  approach  comes  from  conductivity  measurements  of 

structured  regions,  a  topic  that  is  discussed  further  in 

the  next  section. 
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2-3  MICROSTRUCTURE  CONDUCTIVITY  AND  CLOSURE  RELATIONS 

As  seen  in  the  previous  section  the  essential  difference 
between  classical  techniques  and  SCENARIO  microstructure 
techniques  is  contained  within  the  conductivity  relation 
between  the  mean  current  and  the  electric  field.  This 
difference  can  be  seen  by  comparing  equations  (2.12) 

(classical)  and  equation  (2.15)  (SCENARIO).  The  refinements 
provided  by  the  microstructure  approach  are  based  upon  the 
premise  that  fine  scale  structure  alters  the  conductivity 
from  the  value  determined  from  the  mean  plasma  density  in  a 
cell  and  that  by  specifying  additional  statistics  of  the 
plasma  fluctuations  the  conductivity  estimate  can  be  improved. 
This  section  investigates  the  relations  which  prescribe  the 
conductivity  in  terms  of  plasma  density  statistics  and  their 
motivation. 

The  stochastic  potential  equation  (the  stochastic  version 
of  equation  (2.9))  is  a  relation  between  the  deterministic 
neutral  wind,  V,  and  the  stochastic  quantities  of  plasma 
density  N  and  electric  field  potential,  <j> : 

V  •  (NV$)  -  V  •  (N  V  X  i  )  (2.26) 

z 

The  quantities  N,  $  and  V  are  functions  of  spatial 
coordinates  x  and  y.  The  quantities  N  and  <j>  are  also 
inhomogeneous  and  the  potential  solving  step  finds  their 
spatially  varying  statistics  for  the  given  neutral  wind  flow. 

The  numerical  approach  to  the  potential  solving  step  is  to 
spatially  grid  the  problem  to  a  resolution  such  that  the 
spatial  variation  of  the  statistical  descriptors  are  well 
resolved  in  the  sampling  grid.  To  a  good  approximation  the 
stochastic  fluctuations  can  be  treated  as  homogeneous  within 
a  sampling  cell.  By  ascribing  an  effective  conductivity  to 
each  cell,  as  determined  from  solutions  of  the  potential 
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equation  with  homogeneous  processes,  the  spatially  varying 
nature  of  the  problem  cam  be  treated  using  classical  techniques. 
The  stochastic  potential  problem  of  interest  thereby  becomes 
the  determination  of  the  effective  conductivity  of  a  region 
with  homogeneous  statistics.  (The  correlation  properties 
of  the  current  with  certain  functions  of  the  plasma  density 
are  also  of  interest  for  higher  order  moment  convection  but 
little  work  has  been  done  in  this  area.) 

In  many  respects  the  homogeneous  stochastic  potential 
problem  for  structured  plasmas  is  related  to  a  problem  that 
should  be  found  in  other  areas  of  physics  and  engineering, 
namely,  the  problem  of  finding  the  conductivity  of  a 
homogeneous  aggregate  of  materials  with  each  material  having 
a  different  conductivity.  The  problem  of  finding  the 
electrical  conductivity  of  a  slab  of  earth  in  which  copper- 
rich  pebbles  are  found  is  an  example  of  a  similar  problem 
and  there  are  potentially  many  other  problems  that  would  be 
of  interest  in  related  fields.  Yet,  to  the  author's 
knowledge  the  problem  has  not  been  extensively  studied 
which  may  be  somewhat  surprising  since  the  problem  can 
also  be  thought  of  as  the  static  equivalent  to  the  problem 
of  propagation  through  inhomogeneities. 

Solutions  to  the  stochastic  potential  equation  are 
difficult  to  obtain  even  for  the  case  of  homogeneous  statistics. 
There  are  two  special  situations  in  which  solutions  are 
available  and  both  involve  structure  stratification.  These 
cases  are  random  density  profiles  stratified  with  gradients 
parallel  to  the  wind.  Figure  2-3,  and  stratified  with 
gradients  perpendicular  to  the  wind.  Figure  2-4.  Because  of 
the  somewhat  limited  nature  of  the  configurations  which  have 
analytical  solutions  the  application  of  numerical  techniques 


has  proven  useful  for  extending  understanding  of  stochastic 
potential  equation  solutions. 

Consider  the  first  situation  in  which  an  analytical 
solution  exists  as  illustrated  in  Figure  2-3.  This  situation 
could  easily  arise  in  the  simulation  of  smooth  barium  ion 
clouds  that  have  steepened  to  a  point  where  the  actual 
density  gradients  are  unresolvable  by  the  sampling  grid  and 
thus  must  be  treated  as  microstructure.  With  the  wind 
applied  so  that  currents  are  driven  along  the  contours  of 
constant  density,  there  are  no  electric  field  fluctuations 
produced.  The  electric  field  is  constant  and  independent 
of  the  density  fluctuations.  The  average  current  in  the 
cell  is  consequently  given  by 

<J>  -  <N>  (<E>/B  +  V  X  iz)  (2.26) 

and  the  conductivity  is  thus  proportional  to  the  mean  electron 
density.  This  situation  corresponds  to  classical  conductivity 
relations. 

The  other  situation  with  analytical  solution  has  random 
density  fluctuations  with  gradients  perpendicular  to  the  wind 
as  shown  in  Figure  2-4  and  is  similar  to  the  situation  in  late¬ 
time  barium  clouds  in  which  long  sheets  have  developed  parallel 
to  the  wind.  In  this  situation  the  density  is  constant  in  the 
direction  of  the  wind  and  consequently  the  wind  driven  current 
in  the  perpendicular  direction  is  also  constant  in  the  wind 
direction.  Since  the  wind  current  is  divergence  free  it  can 
be  concluded  that  the  current  is  constant.  Since  the  current 
is  constant  it  is  uncorrelated  with  the  fluctuations  in  the 
density.  As  shown  in  Section  2-2  this  fact  leads  to  a  con¬ 
ductivity  that  is  proportional  to  1/<1/N>.  This  form  of 
the  conductivity  is  basically  the  form  used  in  SCENARIO  with 
some  further  approximation. 


31 


N 


<v4>> 


4 


J  - ► 


V 


O  B 


•  Allows  solution  with  $  being  independent  of  N 

•  <J>  =  N  «v<{<>  +  V  xiz) 

•  C  =<N> 

•  This  form  of  conductivity  arises  whenever  N  and  v4> 

ARE  UNCORRELATED 


rigure  2-3.  Structure  stratified  with  gradients 
parallel  to  the  winds. 
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•  Allows  solution  with  J  =  constant.  Thus  J  and  N 

INDEPENDENT 

•  J  -  <l/N>‘1  «v+>  +  V  xiz) 

•  c  =  ci/Nr1 

•  This  form  of  the  conductivity  arises  whenever  J  and 

1/N  ARE  UNCORRELATED 


Figure  2-4.  Structure  stratified  with  gradients 
perpendicular  to  the  wind. 


Since  analytical  solutions  to  the  homogeneous  stochastic 
potential  equation  are  so  limited  numerical  calculational 
results  have  been  used  to  investigate  the  applicability  of 
the  various  conductivity  expressions.  Calculational  results 
provide  the  mean  potential  drop  plus  wind  driver  across  a 
cell  and  the  current  density  at  the  cell  boundary.  The  ratio 
of  these  two  is  the  mean  conductivity. 

Shown  in  Figure  2-5  is  one  geometry  that  has  been 

numerically  investigated.  The  figure  shows  the  upper  half 

of  a  symmetric  problem  of  five  sheet- like  structures.  The 

geometry  is  similar  to  the  stratified  case  investigated 

analytically  but  to  investigate  the  effect  of  the  striation 

tips  upon  the  problem  the  sheets  are  terminated  within 

the  calculational  boundary.  From  the  potential  equation 

solution  for  this  problem  it  is  possible  to  evaluate  the 

conductivity  ratio  on  a  column  by  column  basis  and  to 

compare  the  results  with  the  columnar  statistics  such  as 

<N>  and  1/<1/N> — although  it  should  be  remembered  that 

it  is  the  average  of  these  columnar  ratios  that  should  be 

used  for  the  effective  cell  conductivity.  A  comparison 

of  the  numerically  derived  conductivity  ratio  and  the 

columnar  statistics  is  shown  in  Figure  2-6.  The  figure 

3  2 

shows  the  statistical  ratio  <N>  /<N  >  in  addition  to  the 

quantities  <N>  and  1/<1/N>.  It  can  be  seen  that  on  a 

columnar  basis  the  statistic  1/<1/N>  provides  a  generally 

closer  estimate  to  the  calculated  conductivity  than  the 

3  2 

statistic  of  <N>.  In  addition  the  statistic  <N>  /<N  > 
also  gives  a  respectable  approximation  to  the  numerical 
result.  The  average  of  all  the  columnar  results  is  given 
in  Table  2-1.  The  tabulated  results  show  the  expected 
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Figure  2-6.  Mean  columnar  conductivity  (3)  versus 
position.  Also  shown  are  the  columnar 
mean  density  (0) ,  the  inverse  of  the 
columnar  mean  inverse  density  (1) ,  and 
the  asymtotic  closure  approximation  to 
(1)  used  in  3-dimensional  SCENARIO 
routines  (2). 


Table  2-1.  Comparison  of  Numerical  Conductivity 

Calculational  Results  with  Statistical 


Estimators 

FINGERS 

RODS 

GAUSSIAN 

CLASSICAL  <N> 

1.000 

1.000 

1.000 

SCENARIO  <1/N>~1 

.647 

.682 

.964 

ASYMPTOTIC 

CLOSURE  <N>3/<N2 

>  .640 

.634 

.966 

\ 


NUMERICAL  RESULT  C 


632 


726 


982 


result  that  the  1/<1/N>  SCENARIO  conductivity  provides  a 
better  estimate  to  the  actual  conductivity  than  the  classical 
<N>  value.  Note  that  the  current  pattern  the  magnitude  of 
which  is  shown  in  Fig.  2-7  is  not  smooth  everywhere  in  the 
cell  because  of  the  presence  of  the  striation  tips.  Neverthe¬ 
less,  it  appears  that  it  is  not  crucial  that  the  current  be 
constant  for  the  conductivity  result  to  apply. 

In  the  above  numerical  example  the  structure  is  highly 
stratified  in  the  direction  that  favors  the  SCENARIO  conduc¬ 
tivity  form.  To  remove  the  effects  of  structure  elongation 
from  the  conductivity  comparison  the  conductivity  for  a 
random  distribution  of  nearby  cylindrical  rods  as  shown  in 
Figure  2-8  was  calculated.  The  results  tabulated  in  Table 
2-1  indicate  the  expected  result  that  the  conductivity  does 
not  agree  quite  as  well  as  the  previous  case  with  the 
1/<1/N>  SCENARIO  value.  It  is  surprising  perhaps  that 
they  are  still  much  closer  to  the  SCENARIO  value  than  the 
classical  value.  Thus  use  of  the  SCENARIO  conductivity 
approximation  seems  to  be  justified  for  structured  regions 
in  general  so  long  as  the  density  f lucturations  are  not 
strongly  stratified  with  gradients  parallel  to  the  wind. 

The  numerical  conductivity  examples  shown  in  Figures 
2-5  and  2-7  can  be  closely  described  as  two  dimensional 
shot  noise  processes  (Papoulis,  1965,  p.  288)  in  that  they 
consist  of  a  random  placement  of  similar  shapes.  Numerical 
conductivity  calculations  have  also  been  performed  for  one 
likely  alternative  type  of  stochastic  process,  namely,  the 
two-dimensional  Gaussian  process. 

One  member  of  a  two-dimensional  Gaussian  process 

ensemble  with  isotropic  power  law  power  spectral  density 

is  shown  in  Figure  2-9 .  The  power  law  dependence  used  to 

2  2  -3/2 

generate  the  figure  is  proportional  to  (k  +  k  +  A) 

x  y 
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Figure  2-7.  Contours  of  magnitude  of  current 

density  resulting  from  five  finger 
potential  calculation. 


to  the  outer  scale.  This  power  law  dependence  produces  the 
_2 

k  1-d  in  situ  fluctuation  spectrum  expected  from  many 


structuring  processes. 

The  isotropic  process  may  be  appropriate  for  weakly 
structured  regions  of  plasma.  As  structure  grows  through 
the  gradient  drift  instability  the  non-isotropic  two- 
dimensional  shot  noise  process  probably  becomes  a  better 
descriptor  because  of  the  manner  in  which  structure  tends 
to  form  elongated  finger-like  regions. 

The  process  shown  in  Figure  2-9  is  weakly  structured 
and  has  a  r.m.s.  density  to  mean  density  ratio  of  only  19%. 
The  results  of  the  potential  caclulation  for  this  case  are 
shown  in  Figure  2-10.  Periodic  boundary  conditions  have 
been  used.  The  conductivity  resulting  from  the  calculation 
has  been  compared  with  other  statistical  descriptors  for 
the  conductivity  in  Table  2-1.  The  results  show  that  the 
numerically  calculated  conductivity  is  basically  a  compromise 
between  the  classical  and  SCENARIO  conductivities.  Note 
that  the  asymptotic  closure  conductivity  is  found  to  be 
very  close  to  the  SCENARIO  conductivity.  The  results,  thus, 
show  that  the  SCENARIO  and  asymptotic  closure  conductivity 
approximations  are  no  worse  (and  only  marginally  better)  than 
the  classical  conductivity  for  the  case  of  weak  isotropic 
Gaussian  structure. 

In  all  three  numerical  conductivity  evaluations  the 
3  2 

statistic  <N>  /<N  >  is  seen  to  provide  an  approximation 
that  is  of  comparable  quality  as  the  <1/N>  approximation. 
This  conductivity  relation  has  been  labelled  the  asymptotic 
closure  relation  and  it  has  been  used  extensively  in  the 
SCENARIO  algorithms  to  be  presented.  It  is  instructive 
perhaps  to  discuss  the  origin  of  its  name. 
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The  term  closure  relation  originates  in  past  work 

which  attempted  to  explicitly  calculate  the  statistic  <1/N>. 

This  work  recognized  that  the  statistic  <1/N>  can  be 

convected  from  initial  conditions  and  flow  patterns  in  the 

same  manner  as  <N>  is  convected.  Where  the  convection  of 

<N>  requires  knowledge  of  <1/N>,  the  convection  <1/N> 

2  2 

requires,  analogously,  <1/N  >.  The  statistic  <1/N  >  can 
also  be  convected  provided  <1/N^>  is  known  and  in  general 
the  statistic  <1/N^>  cam  be  convected  provided  that 
<l/N^+^>  is  known,  for  all  p  >  0.  Numerical  treatments 
of  the  problem  require  a  relation  which  terminates  this 
infinite  regression  and  this  relation  was  appropriately 
called  the  closure  relation.  Techniques  were  also  developed 
in  this  early  work  which  did  not  convect  <1/N>  specifically 
but  instead  approximated  it  in  terms  of  other  statistics 
and  the  label  'closure  relation'  was  expanded  to  accommodate 
these  relations.  A  closure  relation  is  thus  any  relation 
between  known  statistics  to  approximate  the  unknown 
statistic  in  the  problem,  including  <1/N>. 

The  statistic  <1/N>  is  directly  evaluatable  from 
the  first  order  statistics  of  N.  Ideally,  one  could  evaluate 
all  the  positive  moments  of  N  and  then  apply  the  statistical 
moment  theorem  to  find  the  first  order  probability  density 
function  and  the  statistic  <1/N>.  In  practice,  however, 

<1/N>  is  approximated  in  one  of  two  ways.  First,  by  assuming 
a  multi-parameter  first  order  density,  such  as  the  "two-level 
variable  epsilon"  densities  used  in  past  work  and  by  then 
using  the  known  statistics  to  determine  the  density  parameters; 
the  first  order  density  and  <1/N>  can  be  found.  Alterna¬ 
tively,  a  relation  for  <1/N>  in  terms  of  the  known 
statistics  can  be  stated  specifically  with  a  family  of  first 
order  density  functions  being  implicitly  assumed.  The 


asymptotic  closure  relation  is  a  closure  relation  of  the 
latter  type  that  was  formulated  such  that  for  both  large  and 
small  amounts  of  structure  the  result  would  'asymptote'  to 
the  value  of  1/<1/N>.  Its  simple  form  in  terms  of  first 
and  second  moments  make  the  relation  attractive  for  numerical 
use. 

The  asymptotic  closure  approximation  was  originally 
developed  as  and  may  be  thought  of  as  an  approximation  to 
1/ < 1/N > .  However,  since  the  statistic  1/<1/N>  is  itself 
only  an  approximation  to  the  conductivity  it  is  conceivably 
more  straightforward  to  regard  the  asymptotic  closure  relation 
as  an  alternative  approximation  to  conductivity.  In  this 
vein  it  is  felt  that  future  work  should  concentrate  more  on 
developing  good  approximations  to  the  conductivity  rather 
than  on  developing  good  closure  approximations  to  <1/N>. 

2-4  NUMERICAL  CONSIDERATIONS 

A  reasonable  question  is  which  of  the  myriad  statistical 
parameters  should  be  used  as  convected  quantities  in  micro¬ 
structure  algorithms.  In  the  past  work,  the  choice  of 
parameters  was  somewhat  ad  hoc.  Recently,  however,  methods 
for  making  choices  between  parameter  sets  have  been  developed. 

In  order  to  perform  a  raicrostructure  calculation  the 
convection  of  at  least  two  statistical  parameters  is  required. 
The  number  of  statistical  quantities  convected  in  a  calculation 
could,  in  principle,  be  quite  large.  It  should  be  kept  in 
mind,  however,  that  the  motivation  for  using  additional 
parameters  must  be  to  improve  conductivity  estimates.  To 
this  end  it  remains  to  be  demonstrated  that  significant 
improvements  in  conductivity  estimates  over  that  provided  by 
a  two-parameter  set  using  the  asymptotic  closure  relation  can 
be  made  by  using  three  or  four  parameter  sets.  Because  of 


the  advantage  of  reduced  storage  requirements  it  is  felt  that 

sets  with  only  two  statistical  parameters  should  be  used  for 

convection  unless  an  improvement  in  conductivity  estimates 

can  be  demonstrated  with  larger  sets.  Thus  one  method  for 

choosing  parameter  sets  is  to  choose  sets  with  only  two  members. 

It  might  be  thought  that  the  set  of  the  mean  density 

2 

<N>  and  the  mean  square  density  <N  >  is  a  reasonable 
parameter  set.  Indeed,  past  work  (Stagat,  et  al. ,  1982)  used 
combinations  of  statistics  including  these  two  for  their 
numerical  algorithms.  Recent  numerical  considerations  have 
revealed  that  it  is  desirable  to  choose  statistical  parameters 
that  are,  in  a  sense,  orthogonal.  It  is  desired  that  the 
values  that  the  parameters  can  take  should  be  independent  of 

2 

each  other.  For  the  <N>  <N  >  set  the  parameter  values  are 

not  independent  of  each  other  because  the  mean  square 

2  2 
<N  >  must  be  greater  than  the  mean  squared  <N>  .  The 

closely  related  set  comprised  of  the  mean,  <N>,  and 
2 

the  variance  £  has  been  found  to  be  more  desirable  for 

2 

use.  The  fact  that  the  variance  E  can  take  values 

independent  of  <N>  and  the  related  fact  that  it  is  less 

susceptible  to  numerical  noise  when  it  is  explicitly 

2 

convected  make  the  <N>  E  set  preferable  to  the  <N> 

<N^>  set. 

The  nature  of  the  numerical  pitfalls  associated  with 
2 

the  <N>  <N  >  set  can  be  illustrated  by  a  classical 

calculation  of  structure  convection.  A  barium  cloud  problem 

2 

has  been  simulated  using  the  <N>  <N  >  pair  and  using 

the  classical  conductivity  proportional  to  <N>.  From  the 

discussion  of  Section  2-2  it  is  known  that  all  the  quantities 
2  2 

<N>  <N  >  and  E  should  convect  with  the  same  velocity. 
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However  explicit  numerical  calculations  performed  in  two 

different  frames  of  reference  bring  out  important  issues. 

The  problem  has  been  initialized  with  a  Gaussian 

spatial  distribution  of  mean  density  and  w^ui  microstructure 

levels  with  r.m.s.  values  of  3%  of  the  mean.  The  initial 

2 

values  of  <N>  and  <N  >  are  shown  in  Figures  2-11  and  2-12. 
With  a  neutral  wind  blowing  in  a  horizontal  direction  toward 
the  right  in  the  figures,  v.he  cloud  has  been  allowed  to  evolve 
300  seconds  in  both  a  frame  at  rest  and  in  a  frame  moving 
at  50  meters  per  second  horizontally  toward  the  right  (1/2 
the  neutral  wind  velocity) .  The  only  source  of  difference 
in  the  two  runs  is  that  all  the  x-component  velocity  values 
used  in  the  moving  frame  convection  step  have  been  reduced 
by  50  meters  per  second  from  the  values  used  in  the  rest 
frame  calculation.  Figures  2-13,  2-14  and  2-15  respectively 
show  the  mean  density,  mean  squre  density  and  the  variance 
for  the  rest  frame  calculation  and  Figures  2-16,  2-17  and 
2-18  respectively  show  the  corresponding  results  from  the 
moving  frame. 

2 

Comparisons  of  the  resulting  values  of  <N>  and  <N  > 

indicate  small  (roughly  5%)  differences  which  by  themselves 

are  not  unreasonable  and  are  thoroughly  acceptable.  In 

2 

contrast,  a  comparison  of  the  resulting  values  of  I  ,  as 

2  2 

determined  from  <N  >  -  <N>  ,  shows  a  striking  disagreement 

in  magnitude  and  morphology.  Further,  neither  of  the  results 

agree  with  the  theoretically  indicated  dependence  which 

specifies  that  the  r.m.s.  fluctuation  level  to  mean  density 

ratio  is  a  constant  3%.  Perhaps  the  most  serious  error  is 

that  the  results  from  both  frames  of  reference  show  regions 
2 

where  Z  is  negative. 
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Figure  2-12.  Initial  density  mean  of  N2.  The 

r.m.s.  density  is  3%  of  the  initial 
density. 


49 


I 


i  i  i  i  i 


i 


x 


Figure  2-13.  Mean  density  after  300  seconds  in 

stationary  frame  of  reference.  The 
plasma  evolution  has  been  simulated 
using  <N>  to  approximate  the 
conductivity. 


2 

14.  Mean  of  N  after  300  seconds 

in  stationary  frame  of  reference. 
This  quantity  has  been  convected 
numerically  with  the  same  velocity 
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Figure  2-15.  The  variance  implied  by  mean  of  H 
minus  the  square  of  the  mean  of  N 
calculated  in  the  stationary  frame 
at  300  seconds.  Note  that  significant 
regions  of  negative  variance  are 
implied. 
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Mean  density  after  300  seconds  in 
translating  frame  of  reference. 

The  frame  of  reference  is  moving 
50  meters  per  second  in  the  direction 
toward  the  100  meter  per  second  wind. 
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Figure  2-18. 


The  variance  implied  by  mean  of  N 
minus  the  square  of  the  mean  of  N 
calculated  in  the  translating  frame 
of  reference.  The  regions  of  negative 
variance  are  morphologically  different 
from  those  of  a  stationary  frame 
indicating  numerical  difficulties  with 
convecting  <N>  and  <N^>. 


The  cause  of  the  erroneous  behavior  of  the  calculated 

variance  can  be  attributed  to  the  fact  that  the  5%  f luctua- 

2 

tions  in  <N>  and  <N  >  become  greatly  exaggerated  because 

of  the  approximate  cancellation  of  the  two  numbers  in  the 
2 

E  calculation.  For  calculations  in  which  structure 

dependent  conductivity  values  are  used  in  place  of  the 

2 

classical  conductivity  the  value  of  E  is  very  important 

in  the  determination  of  the  plasma  flow  and  further  evolution 

2  2 
of  E  .  Convecting  <N>  and  <N  >  explicitly  and  then 

2  2  2 

using  E  *  <N  >  -  <N>  can  produce  even  more  unrealistic 

2 

results  because  the  large  inaccuracies  in  E  resulting 
from  each  time  step  will  feedback  into  the  calculation 
through  the  conductivity  values  and  the  associated  flow 
fields. 

2 

The  major  pitfalls  associated  with  the  <N>  <N  > 

convection  set  are  avoided  if  <N>  and  E^  are  explicitly 

convected.  For  the  classical  conductivity  example  above, 

2 

use  of  the  <N>  E  convection  set  would  produce  acceptable 

2 

results.  The  convection  of  E  is  numerically  no 

2 

different  from  the  convection  of  <N  >  and  consequently 
the  variations  of  the  results  between  the  two  frames  of 
reference  would  have  the  same  5%  level.  This  variation  in 
the  important  E  parameter  is  acceptable  and  argues 
strongly  for  the  use  of  'orthogonal'  statistics  as  convection 
parameters . 

Arguments  for  the  use  of  <;N?  and  as  convection 

parameters  are  that  both  of  these  quantities  are  conserved 
quantities  and  ’that  by  convecting  them  numerical  inaccuracies 
cam  be  guarded  against  by  using  conservation  checks.  The 
variance  parameter  is  in  general  not  a  conserved  quantity  as 
mentioned  in  Section  2.2.  However,  mere  conservation  checks 
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do  not  guarantee  accuracy  as  graphically  pointed  out  by  the 

2 

above  example.  Indeed  use  of  the  <N>  I  convection  set 

has  been  observed  to  come  much  closer  to  conserving  the 
2  2 

quantity  of  <N  >  than  an  <N>  <N  >  calculation  which  adjusts 
2  2 

<N  >  so  that  Z  is  kept  non-negative.  Also  by  not  con¬ 
verting  a  conserved  quantity,  conservation  can  serve  as  a 

numerical  quality  check.  Furthermore,  since  it  is  basically 
2  2 

Z  and  not  <N  >  that  is  of  interest  to  satcom  systems 

.  2 

analysts,  assuring  the  conservation  of  <N  >  is  much  less 

2 

important  than  assuring  the  accuracy  of  Z  . 

2-5  BARIUM  CLOUD  EXAMPLE 

A  numerical  calculation  of  the  two-dimensional  evolution 
of  an  ionospheric  barium  cloud  has  been  performed  using  a 
microstructure  algorithm  implemented  with  the  asymptotic 
closure  conductivity  approximation.  The  problem  geometry 
and  initial  values  are  the  same  as  those  used  for  the 
numerical  example  discussed  in  Section  2-4.  The  barium 
cloud  is  represented  as  an  enhancement  in  the  mean  density 
with  a  Gaussian  spatial  dependence.  The  mean  density  at 
the  enhancement  peak  has  a  value  of  5  and  the  background 
density  has  a  value  of  1.  The  Gaussian  spatial  dependence 
has  an  e-folding  radius  of  10  kilometers.  Microstructure 
is  assumed  to  exist  everywhere  with  an  r.m.s.  level  of  3% 
of  the  mean.  Plots  of  the  initial  density  and  variance  are 
shown  in  Figures  2-19  and  2-20.  The  plot  grid  of  50 
kilometers  by  50  kilometers  is  indicative  of  the  spatial 
region  sampled  by  a  32  by  32  grid.  The  mean  density  and 
the  variance  of  the  density  are  convected  using  the  SCENARIO 
equations  presented  in  Section  2-2.  The  wind  is  assumed  to 
be  horizontal  and  to  the  right  with  a  100  meter  per  second 
speed . 


Distance  (i  100  meters) 


Figure  2-20.  Initial  variance  in  the  barium  cloud 
simulation.  The  standard  deviation 
is  3%  of  the  mean  density. 


The  results  of  the  first  potential  step  calculation  are 
illustrated  in  Figure  2-21.  Because  of  the  free  space 
boundary  conditions  used  in  the  solution  the  grid  boundaries 
do  not  affect  the  motion  of  the  plasma.  The  newly  developed 
free  space  boundary  condition  algorithm  avoids  many  of  the 
deficiencies  associated,  with  previous  potential  solving 
algorithms. 

The  split  step  algorithm  has  been  used  to  advance  the 
density  and  variance  values  through  time.  Essentially  a  ten 
second  interval  has  been  used  for  each  complete  time  step. 

The  predicted  values  at  300'  and  510  seconds  using  a  calcula¬ 
tion  with  a  stationary  frame  of  reference  are  shown  in 
Figures  2-22  thru  2-25 .  The  figures  show  the  physically 
appealing  result  that  the  plasma  structure  region  moves  to 
the  trailing  ( right-most 1  edge  of  the  barium  cloud  and 
grows.  This  behavior  is  in  good  agreement  with  experimental 
observation.  The  growth  rate  indicated  also  seems  to  be 
realistic  if  one  believes  that  barium  cloud  structures  are 
produced  by  gradient  drift  amplification  of  irregularities 
in  the  ambient  and  barium  Pedersen  conductivity  and  that 
these  levels  have  initially  a  1%  r.m.s.  to  mean  density 
ratio.  The  time  at  which  significant-  structure  has  been 
first  observed  in  barium  radio  propagation  data  is  typically 
10-20  minutes.  The  winds  observed  in  these  observations 
produced  an  ion-neutral  slip  velocity  of  roughly  the  same 
magnitude  as  that  observed  in  the  numerical  simulation. 

The  numerical  simulation  data  indicates  a  variance  growth 
rate  of  doubling  every  1 00  seconds.  At  this  rate  the 
structure  variance  becomes  comparable  to  the  square  of  the 
mean  density  level  at  roughly  16  minutes  which  can  be  con¬ 
sidered  to  be  in  good  agreement  with  observations. 
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Figure  2-21.  Results  of  the  first  potential 

solution  step.  Free  space  boundary 
conditions  have  been  implemented. 
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Figure  2-22.  Mean  density  at  300  seconds. 
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The  above  split-step  calculation  has  also  been  performed 
for  a  translating  frame  of  reference.  The  reference  frame 
is  assumed  to  move  at  one  half  of  the  wind  velocity.  Figures 
2-26  thru  2-29  show  the  mean  and  variances  predicted  at  300 
and  510  seconds.  At  300  seconds  the  mean  and  variance  are 
similar  to  those  predicted  in  the  rest  frame  calculations 
although  there  are  small  differences  in  the  mean  and  the  peak 
variance  is  nearly  50%  greater.  By  510  seconds  the  mean  is 
developing  a  slightly  different  appearance  from  the  rest 
frame  calculation.  The  variance  at  510  seconds  is  a  factor 
of  2  greater.  In  view  of  the  exponential  growth  of  structure 
these  differences  are  not  very  significant.  The  differences 
do  however  illustrate  that  microstructure  calculations  are 
more  sensitive  to  numerical  errors  than  classical  calcula¬ 
tions  and  emphasizes  the  importance  of  the  numerical 
considerations  of  Section  2-4.  It  should  be  noted  however 
that  in  spite  of  the  marked  differences  in  the  variance 
values,  the  variance  growth  rates  are  roughly  the  same  as 
those  of  the  rest  frame  calculation  and  can  also  be 
considered  to  be  in  agreement  with  observed  data. 

2-6  TWO-DIMENSIONAL  AXIALLY  SYMMETRIC  BOMB  EVOLUTION 

WITH  RADIAL  HEAVE  WINDS 

The  two-dimensional  techniques  presented  in  Section 
2-2  can  also  be  applied  to  the  case  of  nuclear  bomb  plasma 
evolving  under  the  influence  of  axisymmetric  radially 
directed  heave  winds.  The  effective  driver  of  plasma 
structure  from  high  altitude  nuclear  detonations  is  the 
Pedersen  conductivity  weighted  heave  wind  (see  Section  3-1) . 
For  high  latitude  bursts  this  heave  wind  typically  consists 
of  two  dominant  components.  One  component  is  a  poleward 
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directed  constant  wind.  The  micros tructure  aspects  of  evolu¬ 
tion  under  this  type  of  wind  component  is  different  only  in 
spatial  scale  from  the  barium  cloud  evolution  illustrated  in 
Section  2-5.  The  other  dominant  heave  wind  component  is  an 
axisymmetric  radially  directed  wind.  The  investigation  of 
microstructure  evolution  under  the  influence  of  this  type  of 
wind  is  presented  in  this  section. 

One  interesting  and  useful  facet  of  an  axisymmetric 
radial  wind  which  will  be  capitalized  upon  in  this  section 
is  that  if  the  bomb  plasma  statistics  are  also  axi symmetric 
then  the  mean  electrostatic  potential  is  zero.  This  fact  is 
a  consequence  of  the  property  that  the  wind  driven  component 
of  the  current  is  divergence- free  and  consequently  no  electric 
fields  are  established.  Solutions  of  the  potential  equation 
become  unnecessary  for  following  the  evolution  of  the  plasma. 

Under  the  condition  of  axisymmetry  in  the  wind  and  the 

plasma  statistics  the  two-dimensional  evolution  problem  may 

be  reduced  to  a  one-dimensional  problem.  The  equations  for 

2 

the  one-dimensional  spatial  evolution  of  <N>  and  £ 
through  time  which  are  produced  using  the  SCENARIO  approxima¬ 
tion  can  be  found  from  the  results  of  Section  2-2  to  be: 


3<N> 

3t 


7  •  <N>  F  V 
c 


V  •  (£2  -  F  <N>2)  V 
c 


-  F  V  •  V  <N>2 
C 

where  =  1  -  1/(<N><1/N>) 
c 


(2.27) 


(2.28) 
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where  the  divergence  and  gradient  operations  are  defined 
for  the  radial  coordinate  r  through 

7  *  A  “  x  Tr  {r  V  (2.29) 

V  q  -  |f  (2.30) 

Equations  (2.27)  and  (2.28)  together  with  the  asymptotic 
closure  approximation,  i.e., 

F  a  E2/(<n>2  +  Z2)  (2.31) 

w 

can  be  numerically  time-stepped  to  simulate  plasma  and 
plasma  structure  evolution. 

The  evolution  of  plasma  and  plasma  structure  has  been 
simulated  using  a  generic  wind  model  and  representative 
wind  conditions  as  illustrated  in  Figure  2-30.  The  initial 
plasma  density  is  a  constant  background  value  of  1  plus  a 
Gaussian  of  magnitude  4  which  has  a  100  kilometer  e-fold 
radius.  The  standard  deviation  of  the  plasma  structure 
has  been  set  to  be  3%  of  the  initial  mean  density.  The 
radial  wind  has  a  peak  value  of  1  km/sec  at  the  100  kilometer 
radius  value  with  a  linear  inner  dependence  (from  zero  value) 
and  a  1/r  outer  dependence.  This  wind  profile  is  assumed 
to  be  constant  throughout  the  simulation.  The  spatial 
sampling  rate  used  in  the  simulation  is  2  samples  per 
kilometer. 

The  initial  conditions  have  been  advanced  15  minutes 
in  time  (900  seconds)  using  a  .1  second  time  step.  Figures 
2-31  through  2-34  show  plots  of  the  mean  and  r.m.s.  density 
at  2  1/2  minute  intervals.  Sturcture  growth  is  seen  to 
begin  immediately.  By  300  seconds  the  structure  has  become 
intense  enough  to  significantly  influence  the  motion  of  the 
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Figure  2-32.  Plasma  density  and  structure  at 
300  seconds. 
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Figure  2-36.  Plasma  density  and  structure  at 
900  seconds. 
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mean  plasma  density.  Very  steep  gradients  in  both  density 
and  structure  occur  by  600  seconds.  By  900  seconds  three 
of  these  steep  gradient  fronts  have  become  established. 

At  900  seconds  the  mean  density  has  formed  a  ring 
around  the  interior  peak.  The  exterior  of  the  ring  is  at  a 
radius  of  230  kilometers.  While  this  exterior  radius 
indicates  significant  dispersal  of  the  bomb  plasma  note  that 
the  dispersal  is  significantly  less  than  Kilb's  (Stoeckly 
and  Kilb,  1981)  free  moving  stick  analysis  would  prescribe 
( =  400  kilometers) . 

At  900  seconds  the  plasma  structure  is  intense  both 
at  the  exterior  of  the  mean  density  ring  and  at  the  exterior 
of  the  center  peak  intensity.  The  motion  of  the  mean 
density  with  the  enhancement  of  plasma  structure  at  the 
exterior  of  the  ring  is  physically  suggestive  of  structure 
and  density  being  carried  by  high  density  finger-like  regions 
extending  through  the  ring.  The  structure  in  the  interior 
region  of  the  plasma  is  suggestive  of  the  complementary  low 
density  fingers  extending  into  the  interior  plasma  region.  \ 

The  results  of  this  microstructure  evolution  using 
the  SCENARIO  approximation  is  physically  appealing.  In 
contrast  note  that  if  classical  conductivity  and  correlation 
assumptions  were  used  then  the  plasma  density  and  structure 
would  not  change  from  their  initial  values.  The  structure 
and  density  would  be  frozen  to  their  initial  location  and 
no  striation  dispersal  would  be  predicted. 

Of  numerical  interest  in  this  calculation  are  the 
fronts  of  steep  gradients  in  mean  density  and  structure. 

Special  care  must  be  taken  in  treating  these  gradients 
numerically  to  avoid  Gibbs  phenomena.  Similar  gradients 
would  be  expected  if  the  problem  is  performed  in  two 
dimensions  and  similar  Gibbs  phenomena  may  result  if  care 


is  not  taken.  These  fronts,  thus,  reinforce  the  points 
alluded  to  in  Section  2-4  that  the  SCENARIO  algorithms  are 
somewhat  more  sensitive  than  classical  algorithms  and  that 
care  must  be  taken  to  reduce  sources  of  numerical  noise. 


SECTION  3 


MICROSTRUCTURE  MODELLING  OF  TRANSVERSE  TRANSPORT  FOR 

SCENARIO 

The  development  of  microstructure  transport  techniques 
has  been  directed  toward  the  implementation  of  transverse 
transport  algorithms  in  the  SCENARIO  code.  The  SCENARIO 
code,  in  most  aspects,  is  an  engineering  version  of  the 
MRC  MELT  code.  The  use  of  microstructure  techniques  is 
one  aspect  in  which  SCENARIO  attempts  to  go  beyond  MELT. 

Its  use  is  motived  by  the  need  to  establish  local  structure 
intensities— quantities  not  directly  given  by  MELT — for 
propagation  effect  calculation. 

Like  MELT,  SCENARIO  invokes  the  equi-potential  approxi¬ 
mation  to  reduce  the  three  dimensional  transport  process  to 
two  simpler  processes,  namely,  parallel-to-B  transport  and 
transport  transverse  to  B.  The  parallel  transport  problem 
is  basically  one  of  calculating  the  compressible  flow  of  the 
plasma  in  the  direction  of  the  geomagnetic  field  and  is 
treated  one-dimensionally  with  no  interactions  in  the 
transverse  direction.  The  algorithms  for  parallel  transport 
in  SCENARIO  are  basically  approximations  of  and  simplifica¬ 
tions  to  MELT  algorithms. 

The  problem  of  transverse  transport  in  a  three- 
dimensional  geometry  is  closely  related  to  the  two- 
dimensional  transport  problem  addressed  in  the  previous 
section.  In  this  section  the  microstructure  techniques 
formulated  in  the  previous  section  are  extended  to  the  three- 
dimensional  transverse  transport  problem.  First,  the  general 
theory  is  discussed  in  Section  3-1  and  then  in  Section  3-2 
the  results  of  one  set  of  algorithms  which  have  been 
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implemented  for  SCENARIO  are  presented.  Algorithms  for  alter¬ 
native  assumptions  about  the  three-dimensional  statistical 
characteristics  are  presented  in  Section  3-3. 

3-1  EXTENSION  OF  MICROSTRUCTURE  THEORY  TO  THREE 

DIMENSIONAL  TRANSVERSE  TRANSPORT 

Classical  algorithms  for  three-dimensional  transverse 
transport  numerically  model  the  process  in  a  two-step  fashion 
similar  to  that  discussed  in  Section  2.  In  one  step,  the 
convection  step,  the  continuity  equation  is  used  with  a 
known  two-dimensional  flow  field  to  produce  updated  three- 
dimensional  density  values.  In  the  second  step,  the 
potential  solving  step,  an  updated  two-dimensional  flow 
field  is  calculated  from  the  updated  three-dimensional 
density  values.  Microstructure  influences  the  results  of 
both  steps — a  fact  which  is  ignored  in  classical  techniques. 

In  this  section  the  general  aspects  of  the  extension  of 
microstructure  modelling  techniques  to  this  transport 
problem  are  discussed. 

The  microstructure  three  dimensional  transverse 
transport  split-step  algorithm  is  outlined  in  Figure  3-1. 

It  is  basically  an  extension  of  the  two-dimensional  micro¬ 
structure  split-step  algorithm  seen  in  Figure  2-2  and  has 
four  basic  steps.  First  the  potential  equation  determines 
the  transverse  flow  field  from  specified  spatial  dependences 
of  the  effective  conductivity.  Next,  the  flow  fields  are 
used  to  advance  in  time  the  statistics  of  the  plasma 
density  on  each  transverse  layer.  The  statistics  of 
the  density  are  then  field  line  integrated  to  determine 
the  statistics  of  the  field  line  integrated  Pedersen 
conductivity.  Finally,  the  spatial  dependence  of  the 
effective  conductivity  is  determined  from  the  statistics  of 
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the  field  line  integrated  Pedersen  conductivity.  If  the 
plasma  is  not  convected  to  regions  of  different  neutral 
density  and  if  the  transverse  wind  on  each  layer  is  the 
same  then  the  field  line  integrated  Pedersen  statistics 
are  conserved  and  can  be  convected  using  the  two- 
dimensional  microstructure  technique.  As  noted  by  the 
dotted  line  path  in  the  figure,  the  two-dimensional 
calculation  may  in  some  instances  serve  as  an  alternative 
path  for  determining  these  field  line  integrated 
quantities  and  thereby  allow  consistency  checks  for 
the  convection  algorithm  (Stagat,  et  al.  1982) . 

The  central  features  of  the  three-dimensional 
microstructure  approach  are  basically  the  same  as  those 
used  to  develop  the  two-dimensional  microstructure 
algorithms  of  Section  2.  For  solving  the  potential 
step  fluctuations  in  the  cross-field  conductivity 
are  assumed  to  influence  the  mean  conductivity  through 
some  specified  relation,  typically,  the  asymptotic 
closure  relation.  For  solving  the  convection  step  the 
mean  flux  is  calculated  from  the  mean  cross-field  currents. 
As  in  the  two-dimensional  case,  the  cross-field  currents 
are  assumed  to  be  uncorrelated  with  density  related 
quantities.  These  features  of  the  approach  are  basically 
the  same  as  those  described  in  Section  2.  The  primary 
issues  involved  in  extending  the  theory  to  three 
dimensions  focus  on  how  to  reduce  the  given  three 
dimensional  structure  quantities  to  two-dimensional 
inputs  to  the  potential  equation  and  how  to  construct 
three-dimensional  fluxes  from  its  two-dimensional 
solutions.  Currently  these  issues  are  resolved 
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through  the  use  of  fairly  broad  based  assumptions  based 
upon  physical  intuition. 

With  the  parallels  between  the  three-dimensional 
transverse  transport  problem  and  the  two-dimensional 
evolution  problem  now  established  it  is  appropriate  to 
review  the  derivation  of  the  three-dimensional  relations 
and  to  focus  on  the  two  issues  of  field  integrated  statistics 
and  of  relating  the  local  flux  to  the  plasma  potential. 

As  noted  above  the  issue  of  determining  field  line 
integrated  statistics  from  local  statistics  arises  in 
connection  with  the  solution  of  the  potential  aquation. 

The  potential  equation  arises  from  the  condition  that 
the  total  transverse  current  at  any  field  line  be 
divergenceless  to  assure  that  charge  will  not  accumulate. 

The  field  line  integrated  transverse  current,  J,  is 
given  by 


where  j  is  the  local  transverse  current  and  where  the 
magnetic  field  is  assumed  to  be  z-directed.  Note  that 
lower  case  quantities  will  be  used  throughout  Section  3 
to  denote  three-dimensional  counterparts  to  upper  case 
2-dimensional  quantities  (with  the  exception  of  *  4>  *  which 
denotes  a  2-dimensional  quantity  even  though  it  is  properly 
lower  case) .  Also  note  that  a  constant  magnitude  z- 
directed  magnetic  field  is  assumed  in  order  to  simplify 
the  discussion.  The  actual  SCENARIO  algor ithms  have 
been  developed  for  a  dipole  coordinate  system  and  include 
complicating  coordinate  scale  factors. 


From  Section  2.1  the  local  transverse  current 
can  be  expressed  as  (eq.  2.1). 


j 


en 

yB 


+  V  x  i2) 


(3.2) 


where  n  is  the  local  plasma  density,  v  is  the  local 
wind  velocity,  E  is  the  electric  field,  y  is  the 
ionic  mobility  in  the  neutrals  and  B  is  the  magnitude 
of  the  magnetic  field.  The  assumption  that  the  field 
lines  are  equipotential  implies  that  the  electric  field 
is  a  two-dimensional  quantity.  Performing  the  field 
line  integration  of  this  expression  for  the  local 
current  produces  the  following  relation  for  J 


J  =  N (E/B  +  V  x  i2)  (3.3) 

where  N  and  V  are  the  field  line  integrated  Pedersen 
conductivity  and  field  line  integrated  neutral  wind 
respectively  defined  by 


N  = 

J 

fsg  dz 
|  yB 

(3.4) 

v.j 

f  env  -  /M 

“Ub”  dz/N 

(3.5) 

The  issue  of  determining  the  field  line  integrated 
statistics  from  the  local  statistics  arises  directly  from 
these  two  integrals.  Because  of  their  importance  some 
comments  about  them  are  appropriate.  The  values  of  N 
and  V  result  from  weighted  field  line  integrals  of  the 
local  plasma  density  and  the  wind  velocity.  For  N  the 
weighting  function  of  most  significance  is  the  neutral 
density  which  enters  strongly  in  the  ionic  mobility  term. 
The  neutral  density  weighting  function  varies  rapidly 
with  altitude.  In  the  300  kilometer  altitude  regime  it 
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changes  by  a  factor  of  ten  every  100  kilometers.  For 
V,  as  with  N ,  the  neutral  density  is  also  a  signifi¬ 
cant  weighting  function  but,  unlike  with  N,  the 
integral  for  V  includes  the  local  plasma  density  n 
as  a  weighting  term.  The  field  line  averaged  neutral 
wind,  V  thus  will  be  correlated  with  the  local  density 
statistics. 

The  potential  equation  used  in  classical 
algorithms  results  from  the  assumptions  that  the 
electric  field  is  irrotational  and  that  the  current  J 
is  divergenceless: 


V  •  J  =  V  •  (N (E/B  +  V  X  i  )  ) 

z 

=  7  •  (N(V<p  +  V  x  i  ) )  =  0 

z 


(3 


This  relation  is  the  same  as  the  two-dimensional  relation 
presented  in  Section  2  except  that  here  N  and  V 
denote  field  line  integrated  quantities. 

The  potential  equation  used  for  microstructure 
theory  is  similar.  It  is  derived  from  the  assumption 
that  the  electric  field  is  irrotational  and  that  the 
mean  current,  <J>,  is  divergenceless  yielding: 


V  •  <J>  =  V*  <N(E/B  +  V  x  iz> 

=  V*  <N (V<J>/B  +  V  x  iz)>  =  0 


(3 


As  in  the  two  dimensional  case  approximations  to  the 
mean  conductivity  relation  are  needed  to  reduce  the 
stochastic  potential  equation  to  a  numerical  equation. 
The  requirement  thus  arises  to  define  a  relation  for 
<J>  in  terms  of  the  statistics  of  N,  V<|> ,  and  V  and 
to  then  determine  these  statistics  from  the  local 


statistics  of  the  plasma  density  n.  Determining  the 
statistics  of  field  line  integrated  quantities  from  the 
statistics  of  the  local  plasma  density  requires  knowledge 
about  the  correlation  properties  of  the  plasma  density 
along  field  lines. 

The  range  of  possible  formulations  can  be 
exemplified  by  considering  the  possible  expressions  for 
the  variance  of  N,  denoted  by  Z2  in  terms  of  the 
variance  of  the  local  plasma  density  n  denoted  by  a  . 

For  the  case  of  uncorrelated  structure  between  adjacent 
grid  cells  in  the  z-direction,  the  relation  for  the 
variance  of  the  integrated  Pedersen  conductivity  N  is 

*  /(Sr/  d2  <3-8> 

where  the  integration  is  understood  to  represent  a  summation 
over  the  grid  cells  in  the  z-direction.  For  the  case  of 
structure  which  is  correlated  or  anti-correlated  from  cell 
to  cell  with  correlation  coefficient  of  1  the  relation  would  v>e 


Z 


(3.9) 


where  1^  is  the  region  of  the  fluctuations  positively 
correlated  with  the  field  line  integrated  content  and 
where  is  the  region  of  negatively  correlated 

fluctuations.  Note  that  if  the  value  of  a  in  the 
region  is  assigned  a  negative  value ,  the  above  relation 
can  be  expressed  more  compactly  as 


dz 


Z 


(3.10) 


Giving  a  positive  or  negative  polarity  to  the  value  of  a 

has  been  found  to  be  useful  in  developing  the  completely 

correlated  algorithm  discussed  in  Section  3.2.  Since 
2 

only  cr  is  required  for  propagation  effects  analysis, 
assuming  negative  values  of  a  does  not  alter  any  of 
the  systems  effects  outputs.  Different  types  of  field 
line  correlation  models  can  thus  lead  to  significantly 
different  dependencies  of  the  variance  of  N  on  the 
local  plasma  variance  and  correspondingly  different 
algorithms  as  a  result. 

Assuming  the  issue  of  relating  local  quantities 
to  field  line  integrated  quantities  has  been  resolved, 
solutions  to  the  potential  equation  can  be  obtained 
numerically  producing  values  of  the  mean  potential, 

<V<J»,  and  thereby  precipitating  the  second  issue.  Given 
the  mean  potential,  what  is  the  flux  field  of  the  local 
plasma  structure. 

The  mean  flux  of  the  local  mean  density  and  the 

mean  flux  of  the  mean  square  density  are  formally  given 

2 

by  <n  U>  and  <n  U>  respectively  where  U  is  the 
plasma  velocity  defined  by  V<t>  x  i  .  The  current  issue 
is  basically  how  to  express  these  quantities  in  terms  of 
the  known  statistics  of  n  and  the  result  of  the 
potential  step  <V<1>>  or  equivalently  <U>.  A  first 
step  in  obtaining  these  relations  may  be  taken  as  in 
Section  2  by  investigating  the  expressions  between  the 
plasma  velocity  and  the  current: 
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J  X  i. 


n  J  x  i. 


<n  U>  =  <- 


->  +  <n  V> 


(3.11) 


(3.12) 


n  J  x  i 


<n  U>  =  <- 


->  +  <n2  V> 


(3.13) 


By  invoking  the  assumption  that  the  plasma  density 
related  quantities  1/N,  n/N,  and  n2/N  are  uncorrelated 
with  J ,  evaluating  < J>  in  terms  of  <U> ,  and  then 
substituting;  the  following  mean  flux  relations  are 
obtained. 

<n  U>  =  <^>  <V>  +  <n  V> 

*  <¥ 

<n2  D>  -  <4><U>  T  <V>  ♦  <n2  V> 

The  basic  problem  which  remains  is  how  to  evaluate  the 
quantities 

2 

<^>/  <v>/  V>,  and  <n2  V> 

N  N 

in  terms  of  <n>,  a  ,  <N>,  E,  and  v.  These 

quantities  can  be  evaluated  explicitly  only  if  something 
is  known  about  how  the  local  quantity  n  is  jointly 
distributed  with  the  field  line  integrated  quantity  N. 
Once  again  the  extension  to  the  three-dimensional 
problem  is  focussed  on  its  field  line  correlation 
aspects . 

These  quantities  have  been  evaluated  for  the 
completely  correlated  situation  as  discussed  in  Section 


3-2.  For  the  uncorrelated  situation,  techniques  for 
their  evaluation  are  suggested  in  Section  3-3.  The 
expressions  arrived  at  are  somewhat  more  complex  than 
would  initially  be  guessed.  Neither  does  <^>  equal 
— ^  for  the  completely  correlated  case  nor  does 

<^>  =  <n>  <4>  for  the  uncorrelated  case.  In  both 
cases  further  assumptions  about  the  fluctuation 
statistics  are  required.  These  assumptions  and  the 
resulting  expressions  are  discussed  more  fully  in  the 
next  two  sections. 

3-2  THREE-DIMENSIONAL  TRANSPORT  FOR  COMPLETELY 

CORRELATED  STRUCTURE 

In  order  to  reduce  the  three-dimensional  transverse 
transport  problem  into  a  two-dimensional  problem  suitable 
for  application  of  the  techniques  advanced  in  Section  2 
it  is  necessary  to  know  how  the  structure  is  correlated  in 
the  direction  parallel  to  the  geomagnetic  field.  In  all 
current  approaches  to  this  problem  the  structure  correlation 
properties  are  provided  through  assumption  only.  This 
section  considers  the  case  in  which  the  structure  is  assumed 
to  be  coherently  correlated  along  the  geomagnetic  field 
direction.  Algorithms  for  use  in  SCENARIO  are  developed 
and  selected  results  from  a  SCENARIO  implementation  are 
shown. 

The  coherent  fluctuation  approximation  is  explicitly 
that  there  exists  a  stochastic  form  function  e(x,y)  with 
zero  mean  and  unit  variance  that  describes  the  fluctuations 
of  the  integrated  Pedersen  conductivity  through: 
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N  =  <N>  +  e (x,y)  E  (3.16) 

as  well  as  the  fluctuations  of  the  local  density  through  the 
similar  relations 

n  =  <n>  +  £  (x,y)  a  (3.17) 

Note  that  while  the  r.m.s.  value  of  the  integrated  Pedersen 
conductivity  E  is  defined  in  its  usual  positive  sense, 
the  r.m.s.  values  of  the  local  density  fluctuations  are 
allowed  to  be  either  positive  or  negative.  A  negative 
value  of  a  denotes  that  fact  that  the  local  fluctuations 
n  are  coherently  anti-correlated  with  the  total  electron 
content. 

Coherent  structure  growth  is  produced  when  a 
structured  two  dimensional  velocity  field  is  imposed  upon 
gradients  of  local  density  that  may  or  may  not  be  structured. 
This  structure  growth  is  correlated  with  a  normalized 
correlation  coefficient  of  plus  or  minus  one  according  to 
the  sign  of  an  effective  mixing  velocity  component  in  the 
direction  of  the  local  gradient,  as  will  be  detailed  later 
in  this  section.  The  important  point  to  note  now  is  that 
this  mechanism  for  the  coherent  growth  of  structure  motivates 
the  use  of  the  sign  of  a  as  well  as  its  magnitude.  Further, 
the  important  additional  information  carried  by  the  sign 
of  sigma  imposes  no  added  storage  requirements  to  the 
algorithms.  Its  use  can  be  construed  as  a  more  efficient 
use  of  available  computational  assets. 

The  basic  quantities  in  the  algorithm  are  <n>  and 
a.  The  two  dimensional  values  used  for  input  to  the 
potential  equation,  <N>  and  E,  are  essentially  the 
neutral  density  weighted  sums  of  <n>  and  a  as 
discussed  in  the  previous  section.  Explicity 


<  N> 


(3.18) 


*  ±T n  <n> 

no^L  o 

z“irZnoa  (3-19) 

o 

where  Nq  is  the  neutral  density  normalization  quantity 
defined  by: 

No  *  Y.  no  <3-20) 

where  the  altitude  variations  of  the  magnetic  field  strength 
with  position  along  the  geomagnetic  field  have  been  dropped 
for  simplicity  of  presentation,  and  where  the  sums  are 
assumed  to  be  along  the  geomagnetic  field. 

From  Section  3-1,  equation  (3.6)  ,  the  field  line 
integrated  transverse  current  J  can  be  written  as 

J  -  N  E/B  +  N  V  X  i  (3.21) 

z 

giving 

<J>  *  <N  E/B>  +  <N  V>  x  i  (3.22) 

While  this  expression  is  useful  for  determining  the  mean 
flux,  essentially  <N  E/B>;  for  solving  the  potential 
equation  an  expression  for  the  mean  current  <J>  in  terms 
of  the  mean  electric  field  and  mean  wind  is  needed.  One 
method  to  obtain  this  needed  expression  is  to  divide  the 
equation  for  J  by  N,  to  take  the  expectation  and  to 
then  invoke  the  SCENARIO  approximation  that  is  uncorrelated 
with  1/N  giving: 


<J><1/N>  -  <E/B>  +  <V>  x  i 


(3.23) 


By  invoking  the  asymptotic  closure  approximation  to  the 
effective  conductivity,  <1/N>“1  ,  the  expression  for  the 
mean  current  becomes: 
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<J>  = 


r>  +  <V> 


+  <N>  ‘ 


(.3.24) 


The  potential  equation  which  is  solved  numerically  results 

from  setting  the  divergence  of  this  equation  to  zero.  Note 

2 

that  <N>  and  Z  are  known  from  their  field  line  inte¬ 
grations  and  that  <E/B>  is  to  be  solved  for.  The  neutral 

wind  term  <V>  can  also  be  expressed  in  terms  of  the  local 
wind  v  and  the  local  structure  quantities  through  the 
form  factor  assumption  as  now  discussed. 

The  expression  for  the  field  line  integrated  neutral 
wind  V  from  Section  3-1  is: 


nv  dz 


(3.25) 


Note  that  V  is  a  weighted  integral  of  the  neutral  wind 
where  the  weighting  functions  are  the  l<cal  neutral  density 
nQ  and  the  local  plasma  density  n.  With  the  form  factor 
assumption  V  can  be  written  as 


<n>  v  dz  +  e 


N  (<N>  +  e  Z 
o 


/  V 


v  dz 


(3.26) 


or  for  simplicity 


F  +  e  F 
n _ a 

<N>  +  e  Z 


(3.27) 


where  Fn  and  are  weighted  fluxes  of  plasma  density  and 

structure  which  have  been  defined  through 


?  -In 

n  J 
■/ ", 


<n>  v  dz/N 
o  o 


a  v  dz/N 
o  o 


(3.28) 


(3.29) 
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The  above  expression  for  the  stochastic  variation  of  V 
includes  the  e  form  factor  in  both  the  numerator  and  the 
denominator.  To  find  <V>  long  division  is  used  to 
isolate  e  into  either  the  numerator  or  the  denominator. 
This  technique  is  also  used  frequently  in  the  evaluation 
of  many  flux  related  quantities  as  discussed  later.  For 
evaluating  the  mean  of  V  the  long  division  technique 
produces 


F  -  F  <N>/E 
n  a  ' 

<N>  +  e  1 


(3.30) 


Once  the  form  factor  is  isolated  to  the  denominator  as 
above  an  expectation  produces  <1/N> ,  thus 

F 

<V>  —  -*r-  +  <~>  F  -  <i><N>  F  VZ  (3.31) 

2.  N  n  N  o 

is  the  expression  which  results  and  which  can  be  used 
directly  in  the  potential  equation. 

For  solving  the  convection  equations  for  <n>  and 

2  2 
a  the  mean  fluxes  defined  by  <n  U>  and  <n  U>  are 

required.  As  discussed  in  Section  3-1  these  fluxes  are 

derived  through  relations  for  the  current  with  the  assumption 

2 

that  both  n/N  and  n  /N  are  uncorrelated  with  J  giving 
relations  for  the  mean  fluxes.  Specifically  for  <n  U>  the 
expression  obtained  is 


<n/N>  is  found  from  the  form  factor  assumption  and  the 
factorization  provided  by  one  long  division: 


<£>  -  <n>  teg 
N  <N>  +  £  Z 


f .  _  a  <N>\  „  1 

’  (<n>  -  -j - )  <  — 


,  a 

e  Z  >  Z 


=  (<n>  -  R  a)  <1/N> 


where  R  is  defined  by 

<N>  1 


R  = 


z  “  z  <  1/n> 


<N>  _ 

Z  Fc 


where  Fc  -  1  -  < i/n>  <N> 


(3.33) 


(3.34) 


The  quantity  <n  V>  can  be  found  similarly  to  be 

,  ,  cfF  <N>  » 

<n  V>  =  <|>^<n>  Fn  -  (aFn  +  Fa<n> - 2 - )RJ  (3.35) 

When  these  quantities  are  substituted  into  the  expression 
for  <n  U>  considerable  cancelling  occurs  giving  the 
following  result 


<n  U>  *  <n>  <U>  +  a  W  (3.36) 

where  W  is  the  velocity  vector  defined  by: 

W  =  R(F0/Z  -  <U>)  (3.37) 

This  result  is  physically  appealing  in  that  it  says  the 
flux  of  the  local  plasma  is  determined  by  the  classical 
term  <n>  <U>  plus  a  correction  dependent  upon  the  amount 
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of  local  structure.  The  presence  of  structure  locally 
can  cause  convection  in  the  direction  of  the  mixing  velocity 
W.  The  term  F /E  can  be  recognized  to  be  a  weighted 
field-line  integral  of  the  neutral  wind  when  the  weighting 
function  is  the  local  structure  amplitude.  Thus  W  can 
be  considered  to  be  R  times  an  effective  structured 
wind  force  slip  velocity.  Note  that  with  the  asymptotic 
closure  approximation  R  also  takes  on  a  simple  form, 
namely , 


R 


Z  <N> 
- — 

<N>  +  l 


7 


(3.38) 


which  ranges  from  zero  to  a  maximum  value  of  1/2. 

In  order  to  derive  an  equation  for  the  convection  of 

2  2 
a  the  mean  flux  of  n  is  needed.  This  quantity,  <n  U>, 

can  be  found  in  a  fashion  similar  to  the  above  with  the 

relation 


<n2  U>  =  <~  J  x  i  >  +  <n2  V> 
N  Z 


=  <"  >  <g>  -  -<V>  +  <n2  v> 
<N> 


(3.39) 


In  addition  to  <V>  which  was  evaluated  above,  expressions 
2  2 

for  <n  /N>  and  <n  V>  are  required  and  can  be  found 
using  the  long  division  technique.  The  algebra  is  con¬ 
siderably  more  complicated  but  as  with  <n  U>  considerable 
cancellation  occurs  in  the  derivation  to  produce 

<n2  U>  *  <n>2  <U>  +  2  cnW  -  c2  |  W  +  ^r—2.  (3.40) 


The  convection  equation  for  the  standard  deviation  of 
the  local  density  can  be  expressed  in  terms  of  the  mean  flux 
of  the  density  and  density  squared  by  the  following 
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Vk- 


2  2 
3 a  _  3<n  >  3<n> 

3t  3t  3 t 


(3.41) 


~  V  •  <n  U>  -  i-  V  •  <n2  U> 

3 1  a  2  a 


(3.42) 


Substituting  the  above  relations  for  <n  U>  and  <n  U> 
and  performing  further  cancellation  of  terms  yields  a 
relation  which  is  physically  appealing,  namely 


H  -  -  ^  7  •  <“2o>  -  «  •  7” 


(3.43) 


where 


is  a  velocity  given  by 


q  -  -f  “  r  w 


(3.44) 


If  the  local  density  gradients  are  zero  this  equation 
reverts  to  a  relation  for  the  conservation  of  structure 
variance  moving  with  velocity  Q,  namely 


=  -  V  •  (a2  Q) 


(3.45) 


Thus  local  structure  moves  in  a  conservative  manner  according 
to  the  two-dimensional  flow  field  given  by  Q,  a  sum  of 
effective  structured  wind  flow  field  and  the  mean  plasma 
flow  field  with  each  of  these  velocity  terms  weighted  by 
structure  related  values.  This  flow  is  independent  of  the 
sign  of  a  and  does  not  alter  it. 

If  the  local  density  gradients  are  not  zero  then 
structure  can  grow  according  to  the  W  •  Vn  term.  This 
term  demonstrates  the  growth  expected  to  be  produced  by  a 
structured  flow  field  in  a  smooth  density  gradient.  As 
noted  earlier  in  this  section  the  growth  term  can  be 
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positive  or  negative  according  to  the  relationship  between 
the  local  density  gradient  direction  and  the  velocity  W. 
Positive  structure  growth  denotes  structure  growing  in  a 
manner  correlated  with  fluctuations  in  the  total  electron 
content.  Negative  growth  denotes  growth  toward  anti- 
correlated  structure.  Note  that  the  growth  term  can 
exist  even  if  the  local  structure  o  a  0;  however,  in 
order  to  assure  a  non- zero  value  of  W  structure  must 
exist  somewhere  along  the  field  line,  Z  ?  0,  in  order  for 
structure  to  grow.  As  with  the  two  dimensional  case  an 
initially  smooth,  a  »  0  everywhere,  plasma  will  produce 
no  structure.  All  of  these  features  which  are  evident  in 
the  convection  equations  very  nicely  characterize  the 
features  of  gradient  drift  structure  growth  and  are 
intuitively  appealing  from  a  physical  standpoint. 

The  potential  and  convection  algorithms  have  been  coded 
into  SCENARIO  and  results  from  a  single  burst  Cheyenne  Mt. 
run  are  available.  Figures  3-2  and  3-3  summarize  the  Cartesian 
form  of  these  relations  for  easy  reference.  Shown  in  Figures 
3-4  and  3-5  are  the  me  am  plasma  density  and  the  standard  devia¬ 
tion  of  the  microstructure  used  to  initialize  the  SCENARIO 
plasma  grid.  The  mean  plasma  density  is  determined  by  the  MHD 
initializer  module.  The  standard  deviation  has  been  simply 
chosen  to  be  25%  of  the  mean  density.  Figure  3-6  shows  a  con¬ 
stant  altitude  view  of  the  mean  density  near  the  bottom  of  the 
bomb  density  enhancement.  Of  specific  interest  in  this  plot 
is  the  depletion  of  mean  density  near  the  center  and  an  almost 
ring-like  distribution  of  the  peak  value  of  mean  density.  On 
the  interior  of  the  ring  are  density  gradients  which  would  be 
stable  to  the  heave  winds  in  a  two-dimensional  problem;  however, 
in  the  subsequent  evolution  of  this  3-d  problem  the  gradients 
are  the  source  of  important  structure  growth  that  is  anti¬ 
correlated  with  the  total  electron  content  fluctuations. 
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Given  n,  a,  v  (three-dimensional  quantities)  algorithm 
advances  n,  ct  through  time. 


Define 

-  Field  line  integrated  density 

-  Field  line  integrated  structure 

2  -  fs dz 

-  Field  line  integrated  mean  wind  flux  vector 

F  =  [egZdz 
n  J  11 B 

-  Field  line  integrated  structured  wind  flux  vector 


Structure  Strength  Parameter 


-  Mixing  Velocity  Vector 
W  -  R(Fo/I  -  0) 

where  U  is  result  of  potential  calculation. 

Figure  3-2.  Transverse  transport  algorithm  for  completely 
correlated  structure — definitions  in  terms  of 
elemental  parameters. 
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Figure  3-3.  Transverse  transport  algorithm  for 

completely  correlated  structure — equations. 
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Figure  3-4.  Initial  density  for  SCENARIO  plasma 
grid  calculation  for  Cheyenne  Mt. 
calculation. 
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Figure  3-5. 


Initial  standard  deviation  of  plasma 
density.  The  initial  value  is 
arbitrarily  -set  at  25%  of  the  mean 
density  with  all  fluctuations  assumed 
to  be  coherently  correlated  with  the 
total  electron  content  fluctuations. 
Crosses  denote  grid  cell  centers. 


UMPU 


I*  JUM  At  SCIlt 
lit 

CMTOUIt  .  t*tt»AT 
.%m—T 

- -I - p- 


I ALPHA  1 
TIM  Ml 

■  JtTt*AT  .IISMT  .4ATt*AT  ,lltl*ir 
.  UAC*AI  .  1 1IK*A<  .  l*Tt«AA  .  IAAC*t( 

. - -V .  | - 1 - - 


> 

-^NORTH-SOUTH- 


Figure  3-6.  Initial  density  on  a  constant  altitude 
layer.  The  layer  shown  is  the  third 
layer  of  cells  from  the  bottom.  On  all 
plots  the  bold  and  dotted  contours  are 
respectively  the  maximum  and  minimum 
valued  contours. 


Figures  3-7  and  3-8  show  the  mean  plasma  density  and 
the  standard  deviation  of  the  microstructure  at  450  seconds 
after  250  seconds  of  further  plasma  evolution.  During  this 

Q 

interval  the  peak  density  has  decreased  below  10  electrons 
per  cc  as  a  result  of  the  rapid  expansion  of  plasma  along 
field  lines.  Note  that  a  comparison  with  the  MRC  MELT 
Cheyenne  Mountain  calculations  indicates  that  in  this  version 
of  SCENARIO  the  parallel  transport  algorithm  probably  over¬ 
estimates  this  expansion. 

The  transverse  transport  during  this  interval  on  the 
constant  altitude  layer  shown  in  Figure  3-4  is  illustrated 
in  Figures  3-9  and  3-10.  The  transverse  motion  has  been 
governed  by  a  field  line  averaged  neutral  wind  with  a  1-2 
kilometer  per  second  northward  velocity  (as  referred  to  the 
earth's  surface) .  In  response  to  the  wind  the  structure  has 
moved  across  the  gradients  on  the  north  side  of  the  density 
ring.  The  initial  value  of  a  of  +.25  times  the  mean 
decreased  through  zero  to  a  negative  value  indicating  anti¬ 
correlated  structure  growth.  At  this  time  its  magnitude 
begins  to  approach  the  mean  density. 

Figures  3-11  thru  3-14  show  the  mean  plasma  density 
and  standard  deviation  for  the  results  corresponding  to 
Figures  3-7  thru  3-10  obtained  at  900  seconds.  Continued 
parallel  dynamics  are  dijcernible  in  Figures  3-11  and  3-12. 

Of  interest  is  the  enhanced  transverse  transport  visible  in 
the  constant  altitude  plots  of  Figures  3-13  and  3-14. 
Structure  has  grown  to  a  level  of  order  five  times  the  mean 
density  and  has  moved  toward  the  north.  The  mean  plasma 
density  has  moved  with  it.  The  enhanced  northward  mean 
density  transport  is  even  visible  at  lower  altitudes  in 
the  plot  of  Figure  3-11. 

It  is  evident  from  the  above  results  that  three- 
dimensional  transverse  transport  algorithms  that  produce 
physically  appealing  results  and  that  have  a  sturdy 


106 


KWlt 


II  Jll  It 
*1 

COITMI 


ttati 


■Wlr 


IM.ru*  J 

rim  mi 

.mc**r  .IKMr  .UH+4T  .  MK*«r 
.ltll*«T  »M(*«r  .MKMT  .TT***T  ,HHH> 


•H- 


I 


Figure  3-9.  Plasma  density  on  constant  altitude 
layer  at  450  seconds.  This  layer  is 
the  same  layer  as  shown  in  Figure  3-6 
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Figure  3-10. 


Standard  deviation  of  plasma  density  on 
constant  altitude  layer.  This  layer  is 
the  same  layer  as  shown  in  Figure  3-6. 

Note  the  region  of  negative  values  which 
indicate  anti-correlated  structure  growth. 


Figure  3-11.  Plasma  density  at  900  seconds 


microstructure  theory  foundation  can  be  developed  for 
SCENARIO.  In  this  section  a  set  of  transport  algorithms 
have  been  developed  under  the  assumption  that  structure 
along  field  lines  is  coherently  correlated  or  coherently 
anti-correlated.  Because  they  allow  the  structure  to  be 
coherently  anti-correlated  as  well  as  coherently  correlated, 
the  algorithms  developed  herein  represent  a  significant 
advancement  over  previously  described  three-dimensional 
transverse  transport  for  correlated  structure  (Stagat, 

1982) .  When  the  structure  is  not  allowed  to  grow  in  an 
anti-correlated  fashion  unreliable  results  are  obtained 


especially  in  situations  where  anticorrelated  structure 
growth  would  have  occurred.  The  numerical  results  for  a 
Cheyenne  Mountain  burst  SCENARIO  run  with  these  algorithms 


have  been  seen  to  produce  physically  appealing  results. 

3-3  THREE-DIMENSIONAL  TRANSPORT  ALGORITHMS  FOR  UNCORRELATED 

STRUCTURE 

Nuclear  structure  is  believed  to  be  field  aligned 
over  distances  that  are  large  with  respect  to  their  cross¬ 
field  dimensions.  It  may  be  the  case,  however,  that  over 
the  very  large  SCENARIO  grid  lengths  in  the  direction  of 
the  geomagnetic  field  the  fine  scale  structure  may  signifi¬ 
cantly  decorrelate.  At  least  two  physical  processes  could 
cause  structure  to  lose  field  line  coherence,  namely, 
finite  parallel  conductivity  and  image  growth.  Thus  while 
the  investigation  of  uncorrelated  microstructure  transport 
algorithms  is  not  only  motivated  from  the  standpoint  of 
evaluating  both  extremes  of  field  line  coherence;  it  is 


This  section  develops  the  theory  for  an  algorithm 

for  three-dimensional  microstructure  transport  which  could 

be  implemented  into  SCENARIO.  The  basic  quantities  for 

the  algorithm  are  the  local  mean  density  <n>  and  the  local 

2 

microstructure  variance  a  . 

Two  approaches  are  available  for  calculating  the 

structure  related  conductivity  modifications  necessary  for 

solution  of  the  two  dimensional  potential  problem.  One 

approach  is  to  calculate  the  variance  of  the  integrated 

Pedersen  conductivity,  Z2  and  to  use  this  value  with  the 

asymptotic  closure  approximation  discussed  in  Section  2-4 

2 

to  estimate  <1/N>.  To  calculate  Z  an  approximation 
which  is  appropriate  for  the  incoherently  aligned  structure 
is 

Z2»£<noc)2/N2  (3.46) 

where  N^  is  a  neutral  density  normalization  quantity  with 


and  where  nQ  is  the  local  neutral  density.  Once  Z2  is 

determined  the  asymptotic  closure  approximation  gives  <1/N> 

2  2  3 

as  (Z  +  <N>  )/<N>  .  The  potential  equation  can  then  be 
solved  to  find  <U>. 

As  discussed  in  Sections  3-1  and  3-2  <U>  can  be 

used  to  find  the  mean  flux  of  the  density  and  density  squared 

2 

if  the  five  quantities  <n/N>,  <n  /N>,  <V>,  <n  V>  and 
2 

<n  V>  are  known: 


<n  u>  “  <%>  ~^T7n>V2  +  <n  V> 


.  2  ...  n  .  <U>  -  <V>  _ _ 2 

<„  u>  -  <ir>  -  ■71yN>  +  <n  V- 


(.3.48) 

(3.49) 
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In  order  to  evaluate  the  five  quantities  it  is  necessary 

to  make  assumptions  about  the  joint  statistics  of  n  and 

N.  In  Section  3-2  the  form  function  approximation  allowed 

evaluation  of  these  quantities  for  coherently  aligned 

structure.  For  incoherently  aligned  structure  the  key 

assumption  is  that  the  two  stochastic  quantities  n  and 

N  -  n  are  mutually  independent.  For  this  approximation 

2 

evaluation  of  <n/N>  and  <n  /N>  reduces  to  the  evalua¬ 
tion  of 


n  +  y 


and 


n 


n  +  y 


where  y  is  defined  to  be  N  -  n.  Evaluating  the  expressions 
in  this  format  accounts  for  the  fact  that  even  though  the 
structure  in  any  layer  is  uncorrelated  with  all  other  layers 
it  is  still  partially  correlated  with  the  field  line  integrated 
Pedersen  conductivity.  Note  that  <n/N>  is  not  being  assumed 
to  be  <n>  <1/N>. 

The  next  step  in  the  evaluation  is  to  make  assumptions 
regarding  the  statistics  of  the  local  density.  Of  the 
possible  candidates  of  Gaussian,  log  normal,  two- level  and 
others,  the  two-level  distribution  is  numerically  attractive 
and  reproduces  many  of  the  physical  characteristics  believed 
to  be  present.  The  main  problem  with  the  two  level  approxi¬ 
mation  is  that  for  any  given  pair  of  structure  and  mean 
variance,  the  chosen  levels  must  both  be  positive  to  assure 
physically  realistic  quantities.  This  condition  cannot  be 
met  with  an  equally  weighted  two-level  distribution  and 
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consequently  the  probability  density  at  the  two  levels  must 
be  variably  weighted  according  to  an  assumed  function  of 
the  mean  and  variance.  Note  that  this  condition  also  rules 
out  straightforward  use  of  Gaussian  probability  densities. 
While  the  choice  of  the  two-level  density  weighting  function 
is  somewhat  arbitrary,  it  is  felt  that  the  choice  which 
makes  the  local  probability  densities  satisfy  the  asymptotic 
closure  relation  is  an  appropriate  choice. 

Thus  it  is  assumed  that  on  each  layer  n  is 
distributed  according  to 


p(n) 


'1  -  e;  if  n  =  nT 


£ 

0 


if  n  =  njj 
; otherwise 


(3.50) 


where  nL#  nH»  and  e  are  functions  of  the  mean  value  <n> 
and  variance  a2  such  that 


<i> 

n 


a2  +  <n>2 

- 


(3.51) 


Explicitly  these  relations  are 

+/T 


L 

"H 


<n> 


<n> 


-  e 


1  _  1  J  ^  q1 


a  +  4<n>' 


(3.52) 

(3.53) 

(3.54) 


-I 


Given  these  local  density  functions,  evaluation  of 

2 

the  quantities  <n/(n+y)>  and  <n  /(n+y)>  at  each  layer 
is  straightforward.  For  example 

<"/n*y>-*H  -./^f  +  «L -L «•*» 

where  f(y)  is  the  density  of  y.  The  density  of  y  is 
known  if  the  local  densities  are  assumed  to  be  jointly 
independent.  Explicitly  the  quantity  f (y)  is  a  multiple 
convolution  of  the  two-level  density  function  over  all 
the  other  layers  and  is  thus  comprised  of  2L  -  2  Dirac 
delta  functions  of  unequal  weight  where  L  is  the  number 
of  SCENARIO  cells  along  the  geomagnetic  field.  The  term 
<n/(n+y)>  is  thus  explicitly  the  sum  of  2L  terms  that 

can  be  straightforwardly  evaluated  computationally.  The 

2  2 
values  of  <n  /(n+y)>,  <V>,  <n  V>  and  <n  V>  can  be 

evaluated  similarly. 

2  2 
Given  <n/N>,  <n  /N>,  <V>,  <n  V>,  and  <n  V>  the 

mean  flux  of  the  density  and  the  density  squared  can  be 

determined.  In  a  manner  similar  to  that  used  in  Sections 

2 

2-4  and  3-2  the  n  and  n  continuity  equations  can  be 

2 

manipulated  into  an  equation  for  the  convection  of  a  . 

The  approach  has  as  yet  not  been  fully  developed  for 
SCENARIO  algorithms  for  various  reasons.  The  most 
important  reason  is  that  extensive  computation  time  is 
required  for  evaluating  the  terms  like  <n/(n+y)>  over 
the  multi-layer  multi-level  density  function.  Future 
extension  of  this  approach  most  likely  will  need  to 
consider  approximations  to  f(y)  in  these  evaluations. 

It  is  fair  to  conclude  now,  however,  that  a  technique 
for  calculating  three-dimensional  transverse  transport 
of  uncorrelated  structure  has  been  set  upon  a  firm 
theoretical  footing. 
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It  is  appropriate  at  this  time  to  mention  the  other 
approach  available  for  the  evaluation  of  the  mean  conduc¬ 
tivity  1/<1/N>.  The  probability  density  function  for  N 
can  be  found  in  a  manner  similar  to  that  used  for  finding 
f(y)  in  the  above  evaluation.  That  is,  the  probability 
density  of  N  could  be  expressed  as  the  convolution  of  L 
two-level  local  plasma  density  probability  densities 
yielding  a  function  with  as  many  as  2L  unequally  weighted 
Dirac  delta  functions.  Using  this  probability  density 
function  an  alternative  evaluation  for  <1/N>  for  use 
in  the  potential  equation  is  possible.  Using  this  approach 
in  an  algorithm  might  also  be  considered  more  consistent 
than  using  the  asymptotic  closure  assumption  approach.  It 
should  be  noted  in  the  defense  of  the  asymptotic  closure 
approximation  to  <1/N>  in  the  potential  equation  that  to 
a  large  extent  the  consistency  between  it  and  the  approxima¬ 
tion  to  <n/N>  in  the  convection  equation  is  unimportant  in 
comparison  with  the  other  approximations  involved.  Clearly 
though,  if  <1/N>  is  available  in  the  self-consistent  form 
by  using  the  multi-layer  multi-level  statistics,  its  use 
would  seem  to  be  preferable.  Both  approximations  should  be 
considered  during  any  algorithm  verification. 
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SECTION  4 


CONCLUSIONS 

Numerical  simulations  of  high  altitude  nuclear 
explosions  cannot  be  gridded  fine  enough  to  spatially 
resolve  the  structure  in  the  plasma  density.  Classical 
techniques  ignore  the  presence  of  structure  that  cannot 
be  resolved. 

The  importance  of  the  location  and  intensity  of 
plasma  structure  for  the  evaluation  of  RF  propagation 
effects  motivates  the  application  of  stochastic  theory  to 
the  HANE  evolution  problem.  The  transport  of  plasma  and 
plasma  structure  transverse  to  the  magnetic  field  is  closely 
related  to  the  problem  of  the  evolution  of  a  two-dimensional 
plasma.  In  Section  2  stochastic  theory  is  applied  to  the 
2-d  problem  and  modifications  to  the  classical  2-d  split- 
step  evolution  algorithm  are  developed  which  include  micro¬ 
structure  effects.  In  Section  3  the  2-d  algorithms  are 
extended  to  the  HANE  3-d  transveria  transport  problem  for 
implementation  in  the  SCENARIO  code. 

Numerical  calculations  presented  in  Section  2  demon¬ 
strate  that  structured  plasma  has  a  lower  effective  con¬ 
ductivity  than  unstructured  plasma  with  the  same  mean  density. 
Thus,  classical  algorithms  which  ignore  structure  calculate 
somewhat  erroneous  current  patterns  and  plasma  flow  fields. 
Microstructure  algorithms  attempt  to  correct  this  deficiency 
by  convecting  plasma  density  structure  statistics  as  well 
as  the  mean  plasma  density  and  by  using  structure  statistics 
to  estimate  the  local  conductivity. 

The  difference  between  the  classical  convection 
algorithms  and  the  algorithms  implemented  for  convection 
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in  SCENARIO  can  be  understood  on  a  stochastic  basis  as  which 
quantities  are  assumed  uncorrelated.  In  the  classical 
algorithms  fine  scale  density  fluctuations  are  assumed  to 
be  uncorrelated  with  the  electric  field  fluctuations.  In 
the  SCENARIO  algorithms  the  inverse  of  the  density  fluctua¬ 
tions  are  assumed  to  be  uncorrelated  with  the  fluctuations 
in  the  current.  The  former  tends  to  be  the  case  when 
structure  gradients  are  parallel  with  the  wind  while  the 
latter  applies  both  to  the  case  where  the  gradients  are 
predominantly  perpendicular  to  the  wind  and  to  the  case  of 
isotropic  structure. 

A  split-step  algorithm  which  convects  the  mean  plasma 
density  and  the  r.m.s.  density  fluctuation  level  is  developed 
in  Section  2.  This  algorithm  is  applied  to  the  problem  of  2-d 
evolution  of  a  barium  cloud.  The  results  obtained  can  be  con¬ 
sidered  to  be  in  agreement  with  observed  barium  cloud  structur¬ 
ing  times.  The  algorithm  is  also  applied  to  the  problem  of 
radial  winds  through  an  axisyumetric  2-d  structured  plasma. 

In  Section  2  the  choice  of  which  statistical  parameter 
should  be  used  in  convection  algorithms  is  discussed.  Algo¬ 
rithms  which  use  statistical  parameters  which  are  not  independ¬ 
ent  of  each  other,  e.g. ,  mean  density  and  mean  square  density, 
can  produce  unphysical  results  and  should  be  avoided. 

The  2-d  microstructure  theory  is  extended  to  the  3-d 
transverse  problem  in  Section  3.  A  key  consideration  in 
this  extension  is  the  correlation  of  the  fluctuations  in 
the  direction  parallel  to  the  field.  The  correlation 
properties  are  currently  specified  by  assumption.  Two 
orthogonal  assumptions  have  been  used  to  develop  convection 
algorithms.  One  is  the  assumption  of  completely  coherent 
structure  along  the  magnetic  field.  The  other  is  the  assump¬ 
tion  of  structure  that  is  independent  from  grid  cell  to  grid 
cell  in  the  magnetic  field  direction. 
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The  completely  correlated  algorithms  have  been 
programmed  for  use  in  SCENARIO  and  some  SCENARIO  results 
are  presented.  The  results  demonstrate  the  convection  of 
structure  and  the  enhanced  convection  of  the  mean  density 
associated  with  the  structure.  The  growth  of  coherently 
anti-correlated  structure  is  also  demonstrated.  The 
algorithms  presented  represent  a  significant  advance  over 
previously  used  algorithms  because  they  allow  for  anti¬ 
correlated  structure  growth. 

Microstructure  theory  is  a  powerful  tool  which  has 
not  as  yet  been  fully  exploited.  It  is  perhaps  appropriate 
to  point  out  areas  where  further  application  may  prove  useful. 

Image  growth  of  fine  scale  structure  is  currently  not 
modelled.  Image  growth  is  a  consequence  of  the  fact  that 
transverse  currents  are  carried  by  ions  and  parallel  currents 
are  carried  by  electrons.  This  fact  leads  to  plasma  enhance¬ 
ments  and  depletions  associated  with  structure  convection. 

The  growth  of  these  regions  could  be  modelled  based  on 
calculated  values  of  the  divergence  of  the  local  transverse 
current.  Note  that  image  growth  produces  uncorrelated 
structure. 

A  natural  extension  of  the  algorithms  of  Section  3 
is  the  development  of  an  algorithm  for  a  field  line  correla¬ 
tion  function  with  a  fixed  dependence.  It  is  believed 
currently  that  a  decomposition  of  the  multi-layer  density 
fluctuations  into  a  basis  set  of  uncorrelated  functions 
through  a  Gram-Schmidt  orthonormalization  could  be  used  to 
extend  the  uncorrelated  algorithms  to  this  special  case. 

The  convection  algorithms  would  then  be  used  to  advance  the 
statistics  of  the  basis  functions. 

Given  that  algorithms  for  a  fixed  but  arbitrary 
correlation  function  along  the  field  line  can  be  developed 


a  still  more  ambitious  extension  is  suggested.  By  accounting 
for  the  correlated  nature  of  gradient  drift  growth  and  the 
uncorrelated  nature  of  image  growth,  the  correlation  prop¬ 
erties  along  the  field  could  be  dynamically  altered.  This 
sophisticated  approach  could  then  produce  structure  correla¬ 
tion  along  the  field  as  a  function  of  time  in  contrast  to 
merely  relying  on  assumption.  While  this  type  of  calculation 
might  be  of  interest  on  a  one-of-a-kind  investigation  of 
field  line  correlation;  because  of  its  probable  complexity, 
it  would  most  likely  not  be  suited  for  use  in  engineering 
codes  such  as  SCENARIO.  It  is  mentioned  here  primarily  to 
point  out  the  advances  possible  through  the  application  of 
microstructure  theory  to  the  HANE  simulation  problem. 
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